text:
Loss functions

For the computation of loss functions and their gradients we leverage the functionality already implmented in Flux. All of the loss functions from Flux have been imported:

julia> names(CounterfactualExplanations.Losses)
19-element Vector{Symbol}:
:Losses
:binary_focal_loss
:binarycrossentropy
:crossentropy
:ctc_loss
:dice_coeff_loss
:focal_loss
:hinge_loss
⋮
:logitbinarycrossentropy
:logitcrossentropy
:mae
:mse
:msle
:poisson_loss
:squared_hinge_loss
:tversky_loss

Classification

For most classification tasks the default :logitbinarycrossentropy (binary) and :logitcrossentropy should be sufficient. For both choices the package has been tested and works natively. When using other loss functions, some caution is recommended though:

External loss functions

Some margin-based loss functions like hinge loss do not expect inputs in the domain $\mathcal{Y}=\{0,1\}$, but rather $\mathcal{Y}=\{-1,1\}$. In those case one needs to ensure that the training labels $y$ are encoded accordingly. In order to use distance-based loss functions like mean squared error (MSE) loss needs to be computed with respect to probibilities rather than logits. This is currently not supported and we genenerally recommend not to use distance-based loss functions in the classification setting (more on this below).

Regression

At this point CounterfactualExplanations.jl is designed to be used with classification models, since the overwhelming majority of the existing literature on counterfactual explanations is set in this context. By default margin-based loss functions are used and computed with respect to logits (more on this below). To produce counterfactual explanations for regression problems users currently need to binarize the problem: let t denote some target value for the continuous dependent variable y in the regression context, then we could respecify the dependent variable as  = 0 for all y \< t and  = 1 otherwise. In future work we want to add full support for regression problems.

Methodological background

This is a short tutorial on loss functions and gradients typically involved in counterfactual search. It involves more maths than perhaps some of the other tutorials.

General setup

We begin by restating the general setup for generic counterfactual search. Let t ∈ {0, 1} denote the target label, M the model (classifier) and x′ ∈ ℝ^(D) the vector of counterfactual features (we will assume all features are continuous). Then the differentiable optimization problem in algorithmic recourse is generally of the following form

$$$x\prime = \arg \min_{x\prime} \ell(M(x\prime),t) + \lambda h(x\prime)$$$

where ℓ denotes some loss function targeting the deviation between the target label and the predicted label and h(⋅) acts as a complexity penality generally addressing the realism or cost of the proposed counterfactual.

Loss function ℓ

Different choices for ℓ come to mind, each potentially leading to very different counterfactual outcomes. In practice, ℓ is often implemented with respect to the logits a = w^(T)x rather than the probabilities p(y′=1|x′) = σ(a) predicted by the classifier. We follow this convention here, but as we shall see depeding on the label domain this convention does not work well for every type of loss function. Common choices for ℓ in the literature include margin-based loss function like hinge loss and logit binary crossentropy (or log) loss. Some use distance-based loss such as mean squared error loss (MSE).

Hinge loss

With respect to the logits a = wx hinge loss can be defined as follows

$$$\ell(a,t^*)=(1-a\cdot t^*)_+=\max\{0,1-a\cdot t^*\}$$$

where t^(*) is the target label in { − 1, 1}. Since above we defined t ∈ {0, 1} we need a mapping h : {0, 1} ↦ { − 1, 1}. Specifically, we want to plug in h(t) = t^(*) where h(⋅) is just the following conditional:

\begin{aligned} h(t)&=\begin{cases} -1 && \text{if} && t=0 \\ 1 && \text{if} && t=1 \end{cases} \end{aligned}

Then our loss function as function of t can restated as follows:

$$$\ell(a,t^*)=\ell(a,t)=(1-a\cdot h(t))_+=\max\{0,1-a\cdot h(t)\}$$$

The first-order derivative of hinge loss with respect to the logits a is simply

\begin{aligned} \ell'(a,t)&=\begin{cases} -h(t) && \text{if} && a \cdot h(t)<=1 \\ 0 && \text{otherwise.} \end{cases} \end{aligned}

In the context of counterfactual search the gradient with respect to the feature vector is then:

\begin{aligned} && \nabla_{x\prime} \ell(a,t)&= \begin{cases} -h(t)\mathbf{w} && \text{if} && h(t)\mathbf{w}^Tx\prime<=1 \\ 0 && \text{otherwise.} \end{cases} \end{aligned}

In practice gradients are commonly computed through autodifferentiation. In this tutorial we use the Zygote.jl package which is at the core of Flux.jl, the main deep learning library for Julia.

The side-by-side plot below visualises the loss function and its derivative. The plot further below serves as a simple sanity check to verify that autodifferentiation indeed yields the same result as the closed-form solution for the gradient.

h(t) = ifelse(t==1,1,-1)
hinge(a,t) = max(0,1-a*h(t))
default(size=(500,500))
a = -2:0.05:2
p1 = plot(a, [hinge(a,1) for a=a], title="Loss, t=1", xlab="logits")
p2 = plot(a, [gradient(hinge,a,1)[1] for a=a], title="Gradient, t=1", xlab="logits")
p3 = plot(a, [hinge(a,0) for a=a], title="Loss, t=0", xlab="logits")
p4 = plot(a, [gradient(hinge,a,0)[1] for a=a], title="Gradient, t=0", xlab="logits")
plot(p1, p2, p3, p4, layout = (2, 2), legend = false)
savefig(joinpath(www_path, "loss_grad_hinge.png"))

# Just verifying that the formula for the gradient above indeed yields the same result.
function gradient_man(x,w,t)
𝐠 = ifelse(h(t)*w'x<=1, -h(t)*w, 0)
return 𝐠
end;
plot(a, [gradient_man(a,1,1) for a=a], legend=:bottomright, label="Manual", title="Gradient", xlab="logits")
scatter!(a, [gradient(hinge,a,1)[1] for a=a], label="Autodiff")
savefig(joinpath(www_path, "loss_grad_hinge_test.png"))

Logit binary crossentropy loss

Logit binary crossentropy loss loss (sometimes referred to as log loss) is defined as follows:

\begin{aligned} && \ell(a,t)&=- \left( t \cdot \log(\sigma(a)) + (1-t) \cdot \log (1-\sigma(a)) \right) \\ \end{aligned}

where σ(a) is the logit/sigmoid link function.

Once again for the purpose of counter factual search we are interested in the first-order derivative with respect to our feature vector x′. You can verify that the partial derivative with respect to feature x′_(d) is as follows:

\begin{aligned} && \frac{\partial \ell(a,t)}{\partial x\prime_d}&= (\sigma(a) - t) w_d \\ \end{aligned}

The gradient just corresponds to the stacked vector of partial derivatives:

\begin{aligned} && \nabla_{x\prime} \ell(a,t)&= (\sigma(a) - t) \mathbf{w} \\ \end{aligned}

As before implementation below is done through autodifferentiation. As before the side-by-side plot shows the resulting loss function and its gradient and the plot further below is a simple sanity check.

# sigmoid function:
function 𝛔(a)
trunc = 8.0 # truncation to avoid numerical over/underflow
a = clamp.(a,-trunc,trunc)
p = exp.(a)
p = p ./ (1 .+ p)
return p
end

# Logit binary crossentropy:
logitbinarycrossentropy(a, t) = - (t * log(𝛔(a)) + (1-t) * log(1-𝛔(a)))

p1 = plot(a, [logitbinarycrossentropy(a,1) for a=a], title="Loss, t=1", xlab="logits")
p2 = plot(a, [gradient(logitbinarycrossentropy,a,1)[1] for a=a], title="Gradient, t=1", xlab="logits")
p3 = plot(a, [logitbinarycrossentropy(a,0) for a=a], title="Loss, t=0", xlab="logits")
p4 = plot(a, [gradient(logitbinarycrossentropy,a,0)[1] for a=a], title="Gradient, t=0", xlab="logits")
plot(p1, p2, p3, p4, layout = (2, 2), legend = false)
savefig(joinpath(www_path, "loss_grad_log.png"))

# Just verifying that the formula for the gradient above indeed yields the same result.
function gradient_man(x,w,y)
𝐠 = (𝛔(w'x) - y) .* w
return 𝐠
end;
plot(a, [gradient_man(a,1,1) for a=a], legend=:bottomright, label="Manual", title="Gradient", xlab="logits")
scatter!(a, [gradient(logitbinarycrossentropy,a,1)[1] for a=a], label="Autodiff")
savefig(joinpath(www_path, "loss_grad_log_test.png"))

Mean squared error

Some authors work with distance-based loss functions instead. Since in general we are interested in providing valid recourse, that is counterfactual explanations that indeed lead to the desired label switch, using one of the margin-based loss functions introduced above seems like a more natural choice. Nonetheless, we shall briefly introduce one of the common distance-based loss functions as well.

The mean squared error for counterfactual search implemented with respect to the logits is simply the squared ℓ² norm between the target label and a = w^(T)x:

\begin{aligned} && \ell(a,t)&= ||t-a||^2 \end{aligned}

The gradient with respect to the vector of features is then:

\begin{aligned} && \nabla_{x\prime} \ell(a,t)&= 2(a - t) \mathbf{w} \\ \end{aligned}

As before implementation and visualizations follow below.

mse(a,t) = norm(t - a)^2

NOTE: I hinted above that the convention of taking derivatives with respect to logits can go wrong depending on the loss function we choose. The plot below demonstrates this point: for t = 0 the global minimum of the MSE is of course also at 0. The implication for counterfactual search is that for t = 0 the search stops when w^(T)x′ = 0. But at this point σ(w^(T)x′) = 0.5, in other words we stop right at the decision boundary, but never cross it. We will see an example of this below. Key takeaway: carefully think about the choice of your loss function and DON’T use distance-based loss functions when optimizing with respect to logits.

p1 = plot(a, [mse(a,1) for a=a], title="Loss, t=1", xlab="logits")
p2 = plot(a, [gradient(mse,a,1)[1] for a=a], title="Gradient, t=1", xlab="logits")
p3 = plot(a, [mse(a,0) for a=a], title="Loss, t=0", xlab="logits")
p4 = plot(a, [gradient(mse,a,0)[1] for a=a], title="Gradient, t=0", xlab="logits")
plot(p1, p2, p3, p4, layout = (2, 2), legend = false)
savefig(joinpath(www_path, "loss_grad_mse.png"))

# Just verifying that the formula for the gradient above indeed yields the same result.
function gradient_man(x,w,y)
𝐠 = 2*(w'x - y) .* w
return 𝐠
end;
plot(a, [gradient_man(a,1,1) for a=a], legend=:bottomright, label="Manual", title="Gradient", xlab="logits")
scatter!(a, [gradient(mse,a,1)[1] for a=a], label="Autodiff")
savefig(joinpath(www_path, "loss_grad_mse_test.png"))

Example in 2D

To understand the properties of the different loss functions we will now look at a toy example in 2D. The code below generates some random features and assigns labels based on a fixed vector of coefficients using the sigmoid function.

# Some random data:
using Flux, Random, CounterfactualExplanations.Data
Random.seed!(1234)
N = 25
w = [1.0 1.0]# true coefficients
b = 0
xs, ys = Data.toy_data_linear(N)
X = hcat(xs...)
counterfactual_data = CounterfactualData(X,ys')

The plot below shows the samples coloured by label along with the decision boundary. You can think of this as representing the outcome of some automated decision making system. The highlighted sample was chosen to receive algorithmic recourse in the following: we will search for a counterfactual that leads to a label switch.

using CounterfactualExplanations
using CounterfactualExplanations.Models: LogisticModel
M = LogisticModel(w, [b])

Random.seed!(1234)
x = select_factual(counterfactual_data,rand(1:size(X)[2]))
y = round(probs(M, x)[1])
target = ifelse(y==1.0,0.0,1.0)
γ = 0.75

# Plot with random sample chose for recourse
plt = plot_contour(X',ys,M)
scatter!(plt,[x[1]],[x[2]],ms=10,label="", color=Int(y))
Plots.abline!(plt,-w[2]/w[1],b,color="black",label="",lw=2)
savefig(joinpath(www_path, "loss_example.png"))

Now we instantiate different generators for our different loss functions and different choices of λ. Finally we generate recourse for each of them:

# Generating recourse
Λ = [0.0, 1.0, 5.0] # varying complexity penalties
losses = [:hinge_loss, :logitbinarycrossentropy, :mse]
counterfactuals = []
for loss in losses
for λ in Λ
generator = GenericGenerator(;loss=loss,λ=λ)
t = loss == :hinge_loss ? h(target) : target # mapping for hinge loss
counterfactual = generate_counterfactual(x, t, counterfactual_data, M, generator; γ=γ, T=50)
counterfactuals = vcat(counterfactuals, counterfactual)
end
end

The code below plots the resulting counterfactual paths.

1. Complexity penalty (λ): has the expected effect of penalizing long counterfactual paths: as the distance between x and x′ the penalty exerts more and more pressure on the gradient in the opposite direction ∇ℓ. For large choices of λ valid recourse is not attainable.
2. Confidence threshold (γ): note how for both log loss and hinge loss we overshoot a bit, that is we end up well beyond the decision boundary. This is because above we chose a confidence threshold of γ = 0.75. In the context of recourse this choice matters a lot: we have a longer distance to travel (=higher costs for the individual), but we can be more confident that recourse will remain valid. There is of course an interplay between λ and γ.
3. The choice of the loss function matters: the distance-based MSE does NOT work without further ajustments when optimizing with respect to logits, as discussed above.
# Plotting
k = length(counterfactuals)
function plot_recourse(counterfactual, t)
l = string(counterfactual.generator.loss)
l = l[1:minimum([5,length(l)])]
λ = string(counterfactual.generator.λ)
plt = plot_contour(X',ys,M;colorbar=false,title="Loss: $(l), λ:$(λ)")
plt = plot(plt, size=(floor(√(k)) * 350, ceil(√(k)) * 350))
Plots.abline!(plt,-w[2]/w[1],b,color="black",label="",lw=2)
t = minimum([t, total_steps(counterfactual)])
scatter!(plt, hcat(path(counterfactual)[1:t]...)[1,:], hcat(path(counterfactual)[1:t]...)[2,:], ms=10, color=Int(y), label="")
return plt
end
max_path_length = maximum(map(counterfactual -> total_steps(counterfactual), counterfactuals))
anim = @animate for i in 1:max_path_length
plots = map(counterfactual -> plot_recourse(counterfactual, i), counterfactuals)
plot(plots..., layout = (Int(floor(√(k))), Int(ceil(√(k)))), legend = false, plot_title="Iteration: " * string(i))
end
gif(anim, joinpath(www_path, "loss_paths.gif"), fps=5)


text:
hours calculator

Time calculator that allows you to make hourly calculations and time conversions using the h:min:s (hours, minutes, seconds) format.

Hours calculator, online calculus

Summary :

Time calculator that allows you to make hourly calculations and time conversions using the h:min:s (hours, minutes, seconds) format.

hours_calculator online

Description :

The time calculator allows you to make additions, subtractions and conversions on hours or durations expressed in seconds, minutes, hours.

The duration calculator allows you to add hours or durations in h:min:s format. For example, to add 2 hours with 14 minutes, enter hours_calculator(2:00:00+00:14:00).

Subtracting hours or durations

The duration calculator allows you to subtract hours or durations in the format h:min:s. To subtract 2h 13min 5s from 47min , you have to enter hours_calculator(2:13:05-00:47:00).

Converting times

The time calculator automatically converts time into hours, minutes and seconds. For example, to convert 5h into minutes and seconds, you have to enter hours_calculator(5:00:00).

Syntax :

hours_calculator(h:min:s+h:min:s+...+h:min:s)

Examples :

Calculate online with hours_calculator (hours calculator)

format: plain


text:
Lemma 56.93.2. Let $f : T \to S$ be a morphism of schemes.

1. For $K$ in $D((\mathit{Sch}/T)_{\acute{e}tale})$ we have $(Rf_{big, *}K)|_{S_{\acute{e}tale}} = Rf_{small, *}(K|_{T_{\acute{e}tale}})$ in $D(S_{\acute{e}tale})$.

2. For $K$ in $D((\mathit{Sch}/T)_{\acute{e}tale}, \mathcal{O})$ we have $(Rf_{big, *}K)|_{S_{\acute{e}tale}} = Rf_{small, *}(K|_{T_{\acute{e}tale}})$ in $D(\textit{Mod}(S_{\acute{e}tale}, \mathcal{O}_ S))$.

More generally, let $g : S' \to S$ be an object of $(\mathit{Sch}/S)_{\acute{e}tale}$. Consider the fibre product

$\xymatrix{ T' \ar[r]_{g'} \ar[d]_{f'} & T \ar[d]^ f \\ S' \ar[r]^ g & S }$

Then

1. For $K$ in $D((\mathit{Sch}/T)_{\acute{e}tale})$ we have $i_ g^{-1}(Rf_{big, *}K) = Rf'_{small, *}(i_{g'}^{-1}K)$ in $D(S'_{\acute{e}tale})$.

2. For $K$ in $D((\mathit{Sch}/T)_{\acute{e}tale}, \mathcal{O})$ we have $i_ g^*(Rf_{big, *}K) = Rf'_{small, *}(i_{g'}^*K)$ in $D(\textit{Mod}(S'_{\acute{e}tale}, \mathcal{O}_{S'}))$.

3. For $K$ in $D((\mathit{Sch}/T)_{\acute{e}tale})$ we have $g_{big}^{-1}(Rf_{big, *}K) = Rf'_{big, *}((g'_{big})^{-1}K)$ in $D((\mathit{Sch}/S')_{\acute{e}tale})$.

4. For $K$ in $D((\mathit{Sch}/T)_{\acute{e}tale}, \mathcal{O})$ we have $g_{big}^*(Rf_{big, *}K) = Rf'_{big, *}((g'_{big})^*K)$ in $D(\textit{Mod}(S'_{\acute{e}tale}, \mathcal{O}_{S'}))$.

Proof. Part (1) follows from Lemma 56.93.1 and (56.93.1.1) on choosing a K-injective complex of abelian sheaves representing $K$.

Part (3) follows from Lemma 56.93.1 and Topologies, Lemma 33.4.18 on choosing a K-injective complex of abelian sheaves representing $K$.

Part (5) is Cohomology on Sites, Lemma 21.21.1.

Part (6) is Cohomology on Sites, Lemma 21.21.2.

Part (2) can be proved as follows. Above we have seen that $\pi _ S \circ f_{big} = f_{small} \circ \pi _ T$ as morphisms of ringed sites. Hence we obtain $R\pi _{S, *} \circ Rf_{big, *} = Rf_{small, *} \circ R\pi _{T, *}$ by Cohomology on Sites, Lemma 21.19.2. Since the restriction functors $\pi _{S, *}$ and $\pi _{T, *}$ are exact, we conclude.

Part (4) follows from part (6) and part (2) applied to $f' : T' \to S'$. $\square$

Comment #3219 by David Hansen on

In parts (5) and (6), the functors $Rf'_{small,\ast}$ on the right-hand sides should be replaced by $Rf'_{big,\ast}$.

In your comment you can use Markdown and LaTeX style mathematics (enclose it like $\pi$). A preview option is available if you wish to see how it works out (just click on the eye in the toolbar).


text:
# Viterbi algorithm

The Viterbi algorithm is a dynamic programming algorithm for finding the most likely sequence of hidden states – called the Viterbi path – that results in a sequence of observed events, especially in the context of Markov information sources, and more generally, hidden Markov models. The forward algorithm is a closely related algorithm for computing the probability of a sequence of observed events. These algorithms belong to the realm of information theory.

The algorithm makes a number of assumptions.

• First, both the observed events and hidden events must be in a sequence. This sequence often corresponds to time.
• Second, these two sequences need to be aligned, and an instance of an observed event needs to correspond to exactly one instance of a hidden event.
• Third, computing the most likely hidden sequence up to a certain point t must depend only on the observed event at point t, and the most likely sequence at point t − 1.

These assumptions are all satisfied in a first-order hidden Markov model.

The terms "Viterbi path" and "Viterbi algorithm" are also applied to related dynamic programming algorithms that discover the single most likely explanation for an observation. For example, in statistical parsing a dynamic programming algorithm can be used to discover the single most likely context-free derivation (parse) of a string, which is sometimes called the "Viterbi parse".

The Viterbi algorithm was conceived by Andrew Viterbi in 1967 as a decoding algorithm for convolutional codes over noisy digital communication links. For more details on the history of the development of the algorithm see David Forney's article http://arxiv.org/abs/cs/0504020v2. The algorithm has found universal application in decoding the convolutional codes used in both CDMA and GSM digital cellular, dial-up modems, satellite, deep-space communications, and 802.11 wireless LANs. It is now also commonly used in speech recognition, keyword spotting, computational linguistics, and bioinformatics. For example, in speech-to-text (speech recognition), the acoustic signal is treated as the observed sequence of events, and a string of text is considered to be the "hidden cause" of the acoustic signal. The Viterbi algorithm finds the most likely string of text given the acoustic signal.

## Overview

The assumptions listed above can be elaborated as follows. The Viterbi algorithm operates on a state machine assumption. That is, at any time the system being modeled is in some state. There are a finite number of states. While multiple sequences of states (paths) can lead to a given state, at least one of them is a most likely path to that state, called the "survivor path". This is a fundamental assumption of the algorithm because the algorithm will examine all possible paths leading to a state and only keep the one most likely. This way the algorithm does not have to keep track of all possible paths, only one per state.

A second key assumption is that a transition from a previous state to a new state is marked by an incremental metric, usually a number. This transition is computed from the event. The third key assumption is that the events are cumulative over a path in some sense, usually additive. So the crux of the algorithm is to keep a number for each state. When an event occurs, the algorithm examines moving forward to a new set of states by combining the metric of a possible previous state with the incremental metric of the transition due to the event and chooses the best. The incremental metric associated with an event depends on the transition possibility from the old state to the new state. For example in data communications, it may be possible to only transmit half the symbols from an odd numbered state and the other half from an even numbered state. Additionally, in many cases the state transition graph is not fully connected. A simple example is a car that has 3 states — forward, stop and reverse — and a transition from forward to reverse is not allowed. It must first enter the stop state. After computing the combinations of incremental metric and state metric, only the best survives and all other paths are discarded. There are modifications to the basic algorithm which allow for a forward search in addition to the backwards one described here.

Path history must be stored. In some cases, the search history is complete because the state machine at the encoder starts in a known state and there is sufficient memory to keep all the paths. In other cases, a programmatic solution must be found for limited resources: one example is convolutional encoding, where the decoder must truncate the history at a depth large enough to keep performance to an acceptable level. Although the Viterbi algorithm is very efficient and there are modifications that reduce the computational load, the memory requirements tend to remain constant.

## Algorithm

Suppose we are given a Hidden Markov Model (HMM) with states , initial probabilities of being in state and transition probabilities of transitioning from state to state . Say we observe outputs . The state sequence most likely to have produced the observations is given by the recurrence relations:[1]

Here is the probability of the most probable state sequence responsible for the first observations (we add one because indexing started at 0) that has as its final state. The Viterbi path can be retrieved by saving back pointers which remember which state was used in the second equation. Let be the function that returns the value of used to compute if , or if . Then:

The complexity of this algorithm is .

## Example

Consider two friends, Alice and Bob, who live far apart from each other and who talk together daily over the telephone about what they did that day. Bob is only interested in three activities: walking in the park, shopping, and cleaning his apartment. The choice of what to do is determined exclusively by the weather on a given day. Alice has no definite information about the weather where Bob lives, but she knows general trends. Based on what Bob tells her he did each day, Alice tries to guess what the weather must have been like.

Alice believes that the weather operates as a discrete Markov chain. There are two states, "Rainy" and "Sunny", but she cannot observe them directly, that is, they are hidden from her. On each day, there is a certain chance that Bob will perform one of the following activities, depending on the weather: "walk", "shop", or "clean". Since Bob tells Alice about his activities, those are the observations. The entire system is that of a hidden Markov model (HMM).

Alice knows the general weather trends in the area, and what Bob likes to do on average. In other words, the parameters of the HMM are known. They can be written down in the Python programming language:

states = ('Rainy', 'Sunny')

observations = ('walk', 'shop', 'clean')

start_probability = {'Rainy': 0.6, 'Sunny': 0.4}

transition_probability = {
'Rainy' : {'Rainy': 0.7, 'Sunny': 0.3},
'Sunny' : {'Rainy': 0.4, 'Sunny': 0.6},
}

emission_probability = {
'Rainy' : {'walk': 0.1, 'shop': 0.4, 'clean': 0.5},
'Sunny' : {'walk': 0.6, 'shop': 0.3, 'clean': 0.1},
}

In this piece of code, start_probability represents Alice's belief about which state the HMM is in when Bob first calls her (all she knows is that it tends to be rainy on average). The particular probability distribution used here is not the equilibrium one, which is (given the transition probabilities) approximately {'Rainy': 0.57, 'Sunny': 0.43}. The transition_probability represents the change of the weather in the underlying Markov chain. In this example, there is only a 30% chance that tomorrow will be sunny if today is rainy. The emission_probability represents how likely Bob is to perform a certain activity on each day. If it is rainy, there is a 50% chance that he is cleaning his apartment; if it is sunny, there is a 60% chance that he is outside for a walk.

Alice talks to Bob three days in a row and discovers that on the first day he went for a walk, on the second day he went shopping, and on the third day he cleaned his apartment. Alice has a question: what is the most likely sequence of rainy/sunny days that would explain these observations? This is answered by the Viterbi algorithm.

# Helps visualize the steps of Viterbi.
def print_dptable(V):
print "    ",
for i in range(len(V)): print "%7s" % ("%d" % i),
print

for y in V[0].keys():
print "%.5s: " % y,
for t in range(len(V)):
print "%.7s" % ("%f" % V[t][y]),
print

def viterbi(obs, states, start_p, trans_p, emit_p):
V = []
path =

# Initialize base cases (t == 0)
for y in states:
V[0][y] = start_p[y] * emit_p[y][obs[0]]
path[y] = [y]

# Run Viterbi for t > 0
for t in range(1,len(obs)):
V.append()
newpath =

for y in states:
(prob, state) = max([(V[t-1][y0] * trans_p[y0][y] * emit_p[y][obs[t]], y0) for y0 in states])
V[t][y] = prob
newpath[y] = path[state] + [y]

# Don't need to remember the old paths
path = newpath

print_dptable(V)
(prob, state) = max([(V[len(obs) - 1][y], y) for y in states])
return (prob, path[state])

The function viterbi takes the following arguments: obs is the sequence of observations, e.g. ['walk', 'shop', 'clean']; states is the set of hidden states; start_p is the start probability; trans_p are the transition probabilities; and emit_p are the emission probabilities. For simplicity of code, we assume that the observation sequence obs is non-empty and that trans_p[i][j] and emit_p[i][j] is defined for all states i,j.

In the running example, the forward/Viterbi algorithm is used as follows:

def example():
return viterbi(observations,
states,
start_probability,
transition_probability,
emission_probability)
print example()

This reveals that the observations ['walk', 'shop', 'clean'] were most likely generated by states ['Sunny', 'Rainy', 'Rainy'], with probability 0.01344. In other words, given the observed activities, it was most likely sunny when Bob went for a walk and then it started to rain the next day and kept on raining.

The operation of Viterbi's algorithm can be visualized by means of a trellis diagram. The Viterbi path is essentially the shortest path through this trellis. The trellis for the weather example is shown below; the corresponding Viterbi path is in bold:

When implementing Viterbi's algorithm, it should be noted that many languages use Floating Point arithmetic - as p is small, this may lead to underflow in the results. A common technique to avoid this is to take the logarithm of the probabilities and use it throughout the computation, the same technique used in the Logarithmic Number System. Once the algorithm has terminated, an accurate value can be obtained by performing the appropriate exponentiation.

## Java implementation

import java.util.Hashtable;

public class Viterbi
{
static final String RAINY = "Rainy";
static final String SUNNY = "Sunny";

static final String WALK = "walk";
static final String SHOP = "shop";
static final String CLEAN = "clean";

public static void main(String[] args)
{
String[] states = new String[] {RAINY, SUNNY};

String[] observations = new String[] {WALK, SHOP, CLEAN};

Hashtable<String, Float> start_probability = new Hashtable<String, Float>();
start_probability.put(RAINY, 0.6f);
start_probability.put(SUNNY, 0.4f);

// transition_probability
Hashtable<String, Hashtable<String, Float>> transition_probability =
new Hashtable<String, Hashtable<String, Float>>();
Hashtable<String, Float> t1 = new Hashtable<String, Float>();
t1.put(RAINY, 0.7f);
t1.put(SUNNY, 0.3f);
Hashtable<String, Float> t2 = new Hashtable<String, Float>();
t2.put(RAINY, 0.4f);
t2.put(SUNNY, 0.6f);
transition_probability.put(RAINY, t1);
transition_probability.put(SUNNY, t2);

// emission_probability
Hashtable<String, Hashtable<String, Float>> emission_probability =
new Hashtable<String, Hashtable<String, Float>>();
Hashtable<String, Float> e1 = new Hashtable<String, Float>();
e1.put(WALK, 0.1f);
e1.put(SHOP, 0.4f);
e1.put(CLEAN, 0.5f);
Hashtable<String, Float> e2 = new Hashtable<String, Float>();
e2.put(WALK, 0.6f);
e2.put(SHOP, 0.3f);
e2.put(CLEAN, 0.1f);
emission_probability.put(RAINY, e1);
emission_probability.put(SUNNY, e2);

Object[] ret = forward_viterbi(observations,
states,
start_probability,
transition_probability,
emission_probability);
System.out.println(((Float) ret[0]).floatValue());
System.out.println((String) ret[1]);
System.out.println(((Float) ret[2]).floatValue());
}

public static Object[] forward_viterbi(String[] obs, String[] states,
Hashtable<String, Float> start_p,
Hashtable<String, Hashtable<String, Float>> trans_p,
Hashtable<String, Hashtable<String, Float>> emit_p)
{
Hashtable<String, Object[]> T = new Hashtable<String, Object[]>();
for (String state : states)
T.put(state, new Object[] {start_p.get(state), state, start_p.get(state)});

for (String output : obs)
{
Hashtable<String, Object[]> U = new Hashtable<String, Object[]>();
for (String next_state : states)
{
float total = 0;
String argmax = "";
float valmax = 0;

float prob = 1;
String v_path = "";
float v_prob = 1;

for (String source_state : states)
{
Object[] objs = T.get(source_state);
prob = ((Float) objs[0]).floatValue();
v_path = (String) objs[1];
v_prob = ((Float) objs[2]).floatValue();

float p = emit_p.get(source_state).get(output) *
trans_p.get(source_state).get(next_state);
prob *= p;
v_prob *= p;
total += prob;
if (v_prob > valmax)
{
argmax = v_path + "," + next_state;
valmax = v_prob;
}
}
U.put(next_state, new Object[] {total, argmax, valmax});
}
T = U;
}

float total = 0;
String argmax = "";
float valmax = 0;

float prob;
String v_path;
float v_prob;

for (String state : states)
{
Object[] objs = T.get(state);
prob = ((Float) objs[0]).floatValue();
v_path = (String) objs[1];
v_prob = ((Float) objs[2]).floatValue();
total += prob;
if (v_prob > valmax)
{
argmax = v_path;
valmax = v_prob;
}
}
return new Object[]{total, argmax, valmax};
}
}

## Extensions

With the algorithm called iterative Viterbi decoding one can find the subsequence of an observation that matches best (on average) to a given HMM. Iterative Viterbi decoding works by iteratively invoking a modified Viterbi algorithm, reestimating the score for a filler until convergence.

An alternate algorithm, the Lazy Viterbi algorithm, has been proposed recently.[2] This works by not expanding any nodes until it really needs to, and usually manages to get away with doing a lot less work (in software) than the ordinary Viterbi algorithm for the same result - however, it is not so easy to parallelize in hardware.

The Viterbi algorithm [3] has been extended to operate with a deterministic finite automaton in order to quickly generate the trellis with state transitions pointing back at variable amount of history.

## Notes

1. Xing E, slide 11
2. (PDF). Vehicular Technology Conference. December 2002. pp. 371–375. .
3. Luk, R.W.P.; R.I. Damper (1998). . IEEE Trans. Speech and Audio Processing 6 (3): 217-225. .

## References

• Viterbi AJ (April 1967). . IEEE Transactions on Information Theory 13 (2): 260–269. .  (note: the Viterbi decoding algorithm is described in section IV.)
• Forney GD (March 1973). . Proceedings of the IEEE 61 (3): 268–278. .
• M.S. Ryan and G.R. Nudd. (1993). The Viterbi Algorithm. Technical Report University of Warwick RR-238.
• Shinghal, R. and Godfried T. Toussaint, "Experiments in text recognition with the modified Viterbi algorithm," IEEE Transactions on Pattern Analysis and Machine Intelligence, Vol. PAMI-l, April 1979, pp. 184-193.
• Shinghal, R. and Godfried T. Toussaint, "The sensitivity of the modified Viterbi algorithm to the source statistics," IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. PAMI-2, March 1980, pp. 181-185.
• Rabiner LR (February 1989). . Proceedings of the IEEE 77 (2): 257–286. .  (Describes the forward algorithm and Viterbi algorithm for HMMs).
• Feldman J, Abou-Faycal I, Frigo M (2002). . Vehicular Technology Conference 1: 371–375. .
• .


text:
# What is the easiest method to construct a bezier curve for the half part of a symmetric bezier curve?

Consider the following symmetric (about vertical axis) bezier curve constructed by 4 points A, B, C, and D. Now I want to construct another bezier curve (named L) for just either the left or right part of it.

\documentclass[pstricks,border=24pt]{standalone}
\usepackage{pst-eucl}
\begin{document}

\begin{pspicture}[showgrid=false](-3,-1)(3,3)
\pstGeonode
(3,3){A}
(1,1){B}
(-1,1){C}
(-3,3){D}
\psbezier(A)(B)(C)(D)
\psline[linecolor=red](0,3)(0,-1)
\end{pspicture}

\end{document}

Without solving a polynomial equation first, how can I construct the bezier curve L? Answers with PSTricks (preferred) or TikZ or Asymptote are welcome.

Note that I hate a trial and error method!

• Is it specifically this curve and only the half of it or more general? – percusse Dec 29 '13 at 19:56
• @percusse: \psbezier can only accept 4 points at most so @g.kov solution is enough. No need to expand for higher-order bezier curves consisting of more than 4 points. – kiss my armpit Dec 29 '13 at 20:04
• PDF doesn't support higher order anyway. I meant is this specific to these curves or arbitrary curves? – percusse Dec 29 '13 at 20:09
• @percusse: I don't understand you. The given 4 points can be changed to any numbers as long as its curve is symmetrical about vertical axis. :-) – kiss my armpit Dec 29 '13 at 20:11
• See, we have a constraint of symmetricity? Do you have others that's what I asked. – percusse Dec 29 '13 at 20:21

Here is the Asymptote version, which does not use the general Asymptote commands to split a path (guide), it just performs a simple subdivision of the cubic Bezier segment:

// split.asy:

size(5cm);
pair A,B,C,D;
A=(3,3); B=(1,1); C=(-1,1); D=(-3,3);

pair P,B1,C1,B2,C2;

P=(A+D+3(B+C))/8;

B1 = (A+B)/2;
C1 = ((A+C)/2+B)/2;
C2 = (D+C)/2;
B2 = ((D+B)/2+C)/2;

guide g=A..controls B and C..D;
guide gl=A..controls B1 and C1..P;
guide gr=P..controls B2 and C2..D;
draw(g);
draw(gl,deepgreen+1.6bp+opacity(0.5));
draw(gr,red+1.6bp+opacity(0.5));
draw((0,3)--(0,-1),red);

dot(A--B--C--D--P,UnFill);
dot(B1--C1,deepgreen,UnFill);
dot(B2--C2,red,UnFill);
label("$A$",A,SE);
label("$B$",B,S);
label("$C$",C,S);
label("$D$",D,SW);
label("$P$",P,NE);
label("$B_1$",B1,E,deepgreen);
label("$C_1$",C1,E,deepgreen);
label("$B_2$",B2,W,red);
label("$C_2$",C2,W,red);

To get a standalone split.pdf, run asy -f pdf split.asy.

• I really need the converting formulae. – kiss my armpit Dec 29 '13 at 19:29
• @Stiff Jokes: Btw, the Bezier curve is a very interesting object, you might enjoy studying it in a spare time. – g.kov Dec 29 '13 at 19:43

A translated version for @g.kov's Asymptote answer in PSTricks.

\documentclass[pstricks,border=24pt]{standalone}
\usepackage{pst-eucl}
\begin{document}

\begin{pspicture}[showgrid=true](-3,-1)(3,3)
\pstGeonode
(3,3){A}
(1,1){B}
(-1,1){C}
(-3,3){D}
\psbezier(A)(B)(C)(D)
\psline[linecolor=red](0,3)(0,-1)
\nodexn{.125(A)+.125(D)+.375(B)+.375(C)}{R'}
\nodexn{.5(A)+.5(B)}{P'}
\nodexn{.25(A)+.5(B)+.25(C)}{Q'}
\pstGeonode
(P'){P}
(Q'){Q}
(R'){R}
\psbezier[linecolor=red](A)(P)(Q)(R)
\end{pspicture}

\end{document}


text:
문제1476--Jolly Jumpers

1476: Jolly Jumpers

실행시간 제한: 1 Sec  메모리사용 제한: 64 MB
제출: 722  통과: 179
[제출] [채점기록] [묻고답하기]

문제 설명

A sequence of n > 0 integers is called a jolly jumper if the absolute values of the difference between successive elements take on all the values 1 through n-1. For instance,

1 4 2 3

is a jolly jumper, because the absolutes differences are 3, 2, and 1 respectively. The definition implies that any sequence of a single integer is a jolly jumper. You are to write a program to determine whether or not each of a number of sequences is a jolly jumper.

입력 설명

Each line of input contains an integer n < 3000 followed by n integers representing the sequence.

출력 설명

For each line of input, generate a line of output saying "Jolly" or "Not jolly".

입력 예시 Copy

4 1 4 2 3
5 1 4 2 -1 6

출력 예시 Copy

Jolly
Not jolly


text:
For large structural models, the solution to ${\bf K}_T\Delta {\bf U}={\bf R}$ can be the computational bottleneck during an analysis. Although computing speed and algorithms to solve ${\bf K}_T\Delta {\bf U}={\bf R}$ are very good, you still want to make sure the solution happens as quickly as possible, particularly when inside the double nested time integration/equilibrium solution loop of nonlinear dynamic analysis. Pennies in the change jar add up.

The computational cost of solving simultaneous linear equations is proportional to bN2, where N is the number of equations and b is the bandwidth of ${\bf K}_T$. The bandwidth of a matrix is the largest distance, over all rows of the matrix, from the diagonal entry to the farthest column with a known non-zero entry. You know a priori where the non-zero entries in a stiffness matrix will be based on the element connectivity and the equation numbers assigned to each node.

Generally, linear equation solvers will solve faster when the matrix bandwidth is reduced, i.e., when the off-diagonal non-zeros are close to the diagonal. You can’t change N, but you can change b, which is where equation numbering comes in to play. Reducing the bandwidth also decreases the amount of memory allocated for banded and profile matrix storage schemes.

There are three approaches to equation numbering in OpenSees–Plain, AMD, and RCM. I will demonstrate the numbering approaches using the frame with node tags shown below. There are N=22 simultaneous equations of nodal equilibrium.

The Plain numberer assigns equation numbers in ascending order of node tags input by the user, i.e., start at node 1, number the DOFs, then go to node 2, node 3, etc. The node tags don’t have to be sequential, or all positive–the Plain numberer will go from low to high.

For the example frame, the Plain numberer gives the following stiffness matrix topology based on numbering equations in nodal order 1,2,3,…,9,10. The off-diagonals are spread out quite a bit, with matrix bandwidth b=21.

The AMD numberer uses approximate minimum degree ordering to number the equations. Applying AMD to the frame example, the equations are numbered in the sequence of nodes 1,10,5,6,8,2,9,3,7,4. Compared to the Plain numberer, the off-diagonals are packed tighter with matrix bandwidth b=16.

The RCM numberer numbers equations using the reverse Cuthill-McKee algorithm. For the frame example, the order of nodes is 10,6,5,8,3,2,7,4,9,1. Compared to the AMD numberer, the off-diagonals are packed more tightly around the diagonal and the matrix bandwidth is b=13.

Both AMD and RCM can clean up poor choices of node numbering, either input by a user or generated by an automated meshing routine. Unless you want to control the equation numbers, e.g., for demonstration or instructional purposes, there’s no reason not to use AMD or RCM. Although RCM is the default option for equation numbering in OpenSees, RCM is not categorically better than AMD for all problems.

## 4 thoughts on “Reduce Your Bandwidth”

1. Gustavo Araujo says:

Nice post!

Would you recommend running one step of analysis and check the bandwidth of the assembled matrix to check whether AMD or RCM are the optimal choice or is that going too much into detail?

Liked by 1 person

1. Hola Gustavo,
Gracias! That’s probably not worth the effort, but could be useful for large models where you’re going to do nonlinear response history analysis.
PD

Like

2. Prabakaran Kesavan says:

Dear PD
It is a wonderful article. Clear and informative! Thanks for throwing lights on many nuances of OpenSEES. We hope your post would help us a lot and expect more such post from you in future.

Like

This site uses Akismet to reduce spam. Learn how your comment data is processed.


text:
David
Farruggio

# diff

Differences and Approximate Derivatives

### Syntax and Description

 diff(X) calculates differences between adjacent elements of X along the first array dimension whose size does not equal 1: If X is a vector of length m, then Y = diff(X) returns a vector of length m-1. The elements of Y are the differences between adjacent elements of X. $Y = \left[\begin{array}{cccc} X(2)-X(1) & X(3)-X(2) & ... & X(m)-X(m-1) \end{array}\right]$ If X is a nonvector p-by-m matrix, then Y = diff(X) returns a matrix of size (p-1)-by-m, whose elements are the differences between the rows of X. $Y = \left[ \begin{array}{c} X(2,:)-X(1,:) \\ X(3,:)-X(2,:) \\ \vdots \\ X(p,:)-X(p-1,:) \end{array}\right]$ diff(X,n) calculates the nth difference by applying the diff(X) operator recursively n times. In practice, this means diff(X,2) is the same as diff(diff(X)). diff(X,n,dim) is the nth difference calculated along the dimension specified by dim. The dim input is a positive integer scalar.

Back to functions


text:
# Commenting Regular Expressions

I'm trying to comment regular expressions in JavaScript.

There seems to be many resources on how to remove comments from code using regex, but not actually how to comment regular expressions in JavaScript so they are easier to understand.

Unfortunately, JavaScript doesn't have a verbose mode for regular expression literals like some other langauges do. You may find this interesting, though.

In lieu of any external libraries, your best bet is just to use a normal string and comment that:

var r = new RegExp(
'('      + //start capture
'[0-9]+' + // match digit
')'        //end capture
);
r.test('9'); //true

• Interesting but the obvious limit is that it doesn't let you build regex literals. – Denys Séguret Mar 17 '13 at 16:51
• What can you get with regex literals that you can't get with the "new RegExp()" syntax? – plemarquand Mar 17 '13 at 16:52
• @plemarquand Escaping becomes a pain. When you embed a language in a string, you need to be considerate of its special characters. For example, /s\/\d+/ becomes 's\\/\\d+'. Because of that, you need to be careful when building dynamically. Basically, the same trouble you have on the server-side with the language soup. – Zirak Mar 17 '13 at 16:59

While Javascript doesn't natively support multi-line and commented regular expressions, it's easy enough to construct something that accomplishes the same thing - use a function that takes in a (multi-line, commented) string and returns a regular expression from that string, sans comments and newlines.

The following snippet imitates the behavior of other flavors' x ("extended") flag, which ignores all whitespace characters in a pattern as well as comments, which are denoted with #:

function makeExtendedRegExp(inputPatternStr, flags) {
// Remove everything between the first unescaped # and the end of a line
// and then remove all unescaped whitespace
const cleanedPatternStr = inputPatternStr
.replace(/(^|[^\\])#.*/g, '$1') .replace(/(^|[^\\])\s+/g, '$1');
return new RegExp(cleanedPatternStr, flags);
}

// The following switches the first word with the second word:
const input = 'foo bar baz';
const pattern = makeExtendedRegExp(String.raw
^       # match the beginning of the line
(\w+)   # 1st capture group: match one or more word characters
\s      # match a whitespace character
(\w+)   # 2nd capture group: match one or more word characters
);
console.log(input.replace(pattern, '$2$1'));

Ordinarily, to represent a backslash in a Javascript string, one must double-escape each literal backslash, eg str = 'abc\\def'. But regular expressions often use many backslashes, and the double-escaping can make the pattern much less readable, so when writing a Javascript string with many backslashes it's a good idea to use a String.raw template literal, which allows a single typed backslash to actually represent a literal backslash, without additional escaping.

Just like with the standard x modifier, to match an actual # in the string, just escape it first, eg

foo\#bar     # comments go here

// this function is exactly the same as the one in the first snippet

function makeExtendedRegExp(inputPatternStr, flags) {
// Remove everything between the first unescaped # and the end of a line
// and then remove all unescaped whitespace
const cleanedPatternStr = inputPatternStr
.replace(/(^|[^\\])#.*/g, '$1') .replace(/(^|[^\\])\s+/g, '$1');
return new RegExp(cleanedPatternStr, flags);
}

// The following switches the first word with the second word:
const input = 'foo#bar baz';
const pattern = makeExtendedRegExp(String.raw
^       # match the beginning of the line
(\w+)   # 1st capture group: match one or more word characters
\#      # match a hash character
(\w+)   # 2nd capture group: match one or more word characters
);
console.log(input.replace(pattern, '$2$1'));

Note that to match a literal space character (and not just any whitespace character), while using the x flag in any environment (including the above), you have to escape the space with a \ first, eg:

^(\S+)\ (\S+)   # capture the first two words

If you want to frequently match space characters, this can get a bit tedious and make the pattern harder to read, similar to how double-escaping backslashes isn't very desirable. One possible (non-standard) modification to permit unescaped space characters would be to only strip out spaces at the beginning and end of a line, and spaces before a # comment:

function makeExtendedRegExp(inputPatternStr, flags) {
// Remove the first unescaped #, any preceeding unescaped spaces, and everything that follows
// and then remove leading and trailing whitespace on each line, including linebreaks
const cleanedPatternStr = inputPatternStr
.replace(/(^|[^\\]) *#.*/g, '$1') .replace(/^\s+|\s+$|\n/gm, '');
console.log(cleanedPatternStr);
return new RegExp(cleanedPatternStr, flags);
}

// The following switches the first word with the second word:
const input = 'foo bar baz';
const pattern = makeExtendedRegExp(String.raw
^             # match the beginning of the line
(\w+) (\w+)   # capture the first two words
);
console.log(input.replace(pattern, '$2$1'));

In several other languages (notably Perl), there's the special x flag. When set, the regexp ignores any whitespace and comments inside of it. Sadly, javascript regexps do not support the x flag.

Lacking syntax, the only way to leverage readability is convention. Mine is to add a comment before the tricky regular expression, containing it as if you've had the x flag. Example:

/*
\+?     #optional + sign
(\d*)   #the integeric part
(       #begin decimal portion
\.
\d+  #decimal part
)
*/
var re = /\+?(\d*)(\.\d+)/;

For more complex examples, you can see what I've done with the technique here and here.

• +1 for #read above, I'm not repeating this crap (yes, somebody follows your links). – Denys Séguret Mar 17 '13 at 16:50

I would suggest you to put a regular comment above the line with the regular expression in order to explain it.

You will have much more freedom.

• How would a regular comment above the line help in this case: (?:[a-z0-9!#$%&'*+/=?^_{|}~-]+(?:\.[a-z0-9!#$%&'*+/=?^_{|}~-]+)*|"(?:[\x01-\x08\x0b\x0c\x0e-\x1f\x21\x23-\x5b\x5d-\x7f]|\$\x01-\x09\x0b\x0c\x0e-\x7f])*")@(?:(?:[a-z0-9](?:[a-z0-9-]*[a-z0-9])?\.)+[a-z0-9](?:[a-z0-9-]*[a-z0-9])?|\[(?:(?:25[0-5]|2[0-4][0-9]|[01]?[0-9][0-9]?)\.){3}(?:25[0-5]|2[0-4][0-9]|[01]?[0-9][0-9]?|[a-z0-9-]*[a-z0-9]:(?:[\x01-\x08\x0b\x0c\x0e-\x1f\x21-\x5a\x53-\x7f]|\\[\x01-\x09\x0b\x0c\x0e-\x7f])+)$)\$ This is the regex for email – Benjamin Gruenbaum Mar 17 '13 at 16:39
• You put a multi-line comment above. And you explain it chunk per chunk (the same way as Explosion Pills suggests, but above. It is much more convenient if you need to modify the reg ex.). – Jean Mar 17 '13 at 16:40
• @BenjaminGruenbaum You know your regex for email is probably flawed, right ? – Denys Séguret Mar 17 '13 at 16:40
• @dystroy Hah :) I would never write a RegExp for email validating, I was illustrating that RegularExpressions get really unreadable really fast. This one is from regular-expressions.info/email.html – Benjamin Gruenbaum Mar 17 '13 at 16:42
• +1 for the mail+regex. But to come back to the topic, I have been using above-regex-comments for a long time and I find it really efficient. This leaves space to fully describe what is intended and how it is implemented. – Jean Mar 17 '13 at 16:45


text:
My Math Forum Working out numbers which are in both sequences.
 User Name Remember Me? Password

 Algebra Pre-Algebra and Basic Algebra Math Forum

 April 11th, 2017, 01:15 AM #1 Newbie   Joined: Apr 2016 From: England Posts: 13 Thanks: 2 Working out numbers which are in both sequences. This question's really stumped me, and I can't seem to find any tutorials on how to answer it... The nth term of a sequence is 2n + 1 The nth term of a different sequence is 3n  1 Work out the three numbers that are in both sequences and between 20 and 40 Now the three numbers I managed to get are 25, 31 and 37 (as I assumed it meant the output of both of these sequences needs to be the same as well as being in between 20 and 40). For example 2x18+1= 37 as does 3x12+1= 37. Yet I still feel I'm going the wrong way around this. Many thanks in advance Last edited by skipjack; April 11th, 2017 at 03:05 AM.
April 11th, 2017, 01:45 AM   #2
Senior Member

Joined: Sep 2015
From: USA

Posts: 1,791
Thanks: 923

Quote:
 Originally Posted by WarmCustard This question's really stumped me, and I can't seem to find any tutorials on how to answer it... The nth term of a sequence is 2n + 1 The nth term of a different sequence is 3n  1 Work out the three numbers that are in both sequences and between 20 and 40 Now the three numbers I managed to get are 25, 31 and 37 (as I assumed it meant the output of both of these sequences needs to be the same as well as being in between 20 and 40). For example 2x18+1= 37 as does 3x12+1= 37. Yet I still feel I'm going the wrong way around this. Many thanks in advance
$2n+1$ will be all the odd numbers

$3n+1$ will be odd for $n$ even

the numbers with even index in the 2nd sequence between 20 and 40 are

25,31,37 as you discovered.

Last edited by skipjack; April 11th, 2017 at 03:06 AM.

 Tags numbers, sequences, working

 Thread Tools Display Modes Linear Mode

 Similar Threads Thread Thread Starter Forum Replies Last Post Martin Hopf Number Theory 17 April 4th, 2015 04:13 PM wilhelm Number Theory 5 November 30th, 2012 09:23 AM doctorbleachers Algebra 1 May 22nd, 2012 12:35 PM momo Number Theory 19 November 22nd, 2007 09:28 AM wilhelm Applied Math 0 December 31st, 1969 04:00 PM

 Contact - Home - Forums - Cryptocurrency Forum - Top

Copyright © 2018 My Math Forum. All rights reserved.


text:
how would i intergrate these.. how do i know when to use the chain rule when i m integrating these??

$\int xe^{-x^2}$

$\int e^{\frac{-y}{2}}$

$\int e^{2x}$

2. Chain rule for integration? Do you mean integration by parts?

u substitution will be sufficient for all of these.

Originally Posted by Legendsn3verdie
$\int xe^{x^2}$
$u = x^2$
Originally Posted by Legendsn3verdie
$\int e^{\frac{-y}{2}}$
$u=\frac{-y}{2}$
Originally Posted by Legendsn3verdie
$\int e^{2x}$
$u=2x$

Also, don't forget to put dx !

3. Originally Posted by Legendsn3verdie
how would i intergrate these.. how do i know when to use the chain rule when i m integrating these??

$\int xe^{-x^2}~{\color{red}dx}$
substitution, $u = x^2$

$\int e^{\frac{-y}{2}}~{\color{red}dy}$
substitution, $u = - \frac y2$

or use the rule $\int e^{kx}~dx = \frac 1ke^{kx} + C$ for $k \ne 0$

$\int e^{2x}~{\color{red}dx}$
see above

EDIT: Well, gee, thanks a lot Wingless!


text:
# Thom class in homology for defining orientability of vector bundles

A rank $$r$$ real vector bundle $$p : E \to B$$ is said to be orientable if there is a Thom class $$\tau \in H^r (D(E),\partial D(E) ; \mathbb{F})$$ (where $$D(E)$$ is a unit disk bundle and $$\mathbb{F}$$ a field) if its restriction to each fibre $$(D(E)_b,\partial D(E)_b)$$ for $$b \in B$$ is a generator of the group $$H^r ( D(E)_b, \partial D(E)_b ; \mathbb{F})$$.

If we swap the cohomology groups for homology groups (by Poincare duality) does that given an equivalent definition?

Thank you.

• So a generator of $H_n(D(E))$ that is in the image of every induced map $H_n(\partial D(E)_b)\rightarrow H_n(D(E))$, where $b\in B$ and $n=\dim(B)$? This is equivalent, but more difficult to work with. Also it only applies in the special case that $p:E\rightarrow B$ is a smooth Riemannian vector bundle, and even here it is generally easier to work with smooth forms than possibly non-representable homology classes. – Tyrone May 21 at 10:45
• @Tyrone thanks, can you please give me some more detail or a reference about to get from my statement to your statement. Is this Lefschetz duality? I am ok with this working only on a smooth riemannian vector bundle because this is actually what I have. – Nasos Evangelou May 22 at 12:54
• @Tyrone also, the Thom isomorphism in cohomology is $H^k(X,\mathbb{F}) \simeq H^{k+r} (D(E), \partial D(E) ; \mathbb{F})$. Is the correct homology version just $H_k(X,\mathbb{F}) \simeq H_{k+r} (D(E), \partial D(E);\mathbb{F})$ ? Also, can you please tell me where the riemannian structure is needed in what you said? – Nasos Evangelou May 22 at 13:26
• I have just applied Poincare duality, which is induced by taking the cap product with the fundamental class. $D(E)$ is an $(n+r)$-dimensional manifold with boundary, and $B$ is $n$-dimensional. The Thom isomorphism on cohomology is induced by taking the cup product with the Thom class $\tau$. Once you understand how cap and cup products interact you see that the Thom isomorphism in homology is $H_{k+r}(D(E),S(E))\xrightarrow{\cong}H_k(B)$, $x\mapsto \tau\cap x$. – Tyrone May 23 at 9:59
• The maps $H^r(D(E),S(E))\rightarrow H^r(D(E)_b,S(E)_b)$ and $H_{n}(D(E)_b)\rightarrow H_n(D(E))$ are Poincare dual. You don't actually need the Riemannian structure to have a Thom class, but if you are specifically working with the disc and sphere bundles then these need to be defined with a metric. – Tyrone May 23 at 10:27


text:
Power of Two

07/16/2016 Math Bit Manipulation

Question

Given an integer, write a function to determine if it is a power of two.

Solution

Result: Accepted Time: 4 ms

Here should be some explanations.

bool  isPowerOfTwo(const int n) {
return  n > 0 &&n==((-n)&(n));
}

Complexity Analytics

• Time Complexity: $O(1)$
• Space Complexity: $O(1)$


text:
Why do integrals always have a dx?

There's so much confusion around dx, especially among Calc 1 and Calc 2 students. And with good reason. dx doesn't start to feel relevant until later on, when you start working with multivariable functions, multiple integrals, and differential equations.

But that's no reason not to at least address dx in an introductory way. In this video, we'll talk about a few of the simplest explanations for what dx does, what it tells us to do, and what it represents.

When you're just starting out with integrals, for the most part, it's okay to think about dx as just notation. It just comes with the integral symbol, and you don't need to know exactly why it's there in order to know that you're supposed to integrate.

Once you start changing variables, it's important to know that the dx always tells you which variable to integrate. That's why you'll see dx with a function in terms of x, dt with a function in terms of t, and dy with a function in terms of y. You'll also need to understand that your limits of integration (if you're dealing with a definite integral), will need to match the variable in your function and the variable in your dx/dt/dy/etc.

When you get to multivariable calculus and multiple integrals, you'll need to realize that seeing dx dy at the end of the integral tells you to integrate first with respect to x, and then with respect to y.

But most importantly, you'll want to understand that dx represents the differential, or the difference between two values of x. It's the distance between two values of x. So, when you use a Riemann sum, or trapezoidal rule to approximate the area under a curve, delta x is the width of each rectangle (or trapezoid).

When you switch from these approximation methods, and start using the integral instead to find exact area, that delta x has to change to dx in order to match the integral. Essentially what you're doing is saying that delta x represents a larger width, which is why you can only get an approximation. But dx represents and infinitely small distance, which is why it gets paired with the integral, and allows you to find exact area, instead of just an approximation.

And that's why you always need a dx whenever you're using an integral.


text:
It is currently 27 Jun 2017, 22:46

### GMAT Club Daily Prep

#### Thank you for using the timer - this advanced tool can estimate your performance and suggest more practice questions. We have subscribed you to Daily Prep Questions via email.

Customized
for You

we will pick new questions that match your level based on your Timer History

Track

every week, we’ll send you an estimated GMAT score based on your performance

Practice
Pays

we will pick new questions that match your level based on your Timer History

# Events & Promotions

###### Events & Promotions in June
Open Detailed Calendar

# Car Dealership X, the total profit from sales increased by 10 percent

 new topic post reply Question banks Downloads My Bookmarks Reviews Important topics
Author Message
TAGS:

### Hide Tags

Senior Manager
Joined: 25 Feb 2010
Posts: 461
Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

23 Sep 2010, 14:47
2
This post was
BOOKMARKED
00:00

Difficulty:

45% (medium)

Question Stats:

66% (02:29) correct 34% (01:49) wrong based on 88 sessions

### HideShow timer Statistics

Car Dealership X, the total profit from sales increased by 10 percent over the previous years, while the number of cars sold decreased by 10 percent over the previous years. Approximately what was the average percent increase in profit per car over the previous year.

A. 18%
B. 20%
C. 22%
D. 23%
E. 25%
[Reveal] Spoiler: OA

_________________

GGG (Gym / GMAT / Girl) -- Be Serious

Its your duty to post OA afterwards; some one must be waiting for that...

Retired Moderator
Joined: 02 Sep 2010
Posts: 803
Location: London
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

23 Sep 2010, 14:52
1
KUDOS
onedayill wrote:
Car Dealership X, the total profit from sales increased by 10 percent over the previous years, while the number of cars sold decreased by 10 percent over the previous years. Approximately what was the average percent increase in profit per car over the previous year.

18%
20%
22%
23%
25%

Let Profit be x and cars sold be y

Originally Profit/Car = x/y

Now it is 1.1x / 0.9 y = 11/9 (x/y)

Increase in profit per car = ((2/9) (x/y)) / (x/y) * 100 = 22.22%

_________________
Math Expert
Joined: 02 Sep 2009
Posts: 39725
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

23 Sep 2010, 14:59
onedayill wrote:
Car Dealership X, the total profit from sales increased by 10 percent over the previous years, while the number of cars sold decreased by 10 percent over the previous years. Approximately what was the average percent increase in profit per car over the previous year.

18%
20%
22%
23%
25%

Let the initial profit be $$p$$ and initial number of cars sold be $$n$$ --> initial profit per car would be $$\frac{p}{n}$$;

After 1 year profit per car would be $$\frac{1.1p}{0.9n}\approx{1.22*\frac{p}{n}}$$;

General formula for percent increase or decrease, (percent change):

$$Percent \ change=\frac{Change}{Original}*100=\frac{1.22*\frac{p}{n}-\frac{p}{n}}{\frac{p}{n}}*100=22$$.

_________________
Senior Manager
Joined: 25 Feb 2010
Posts: 461
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

23 Sep 2010, 15:00
shrouded1 wrote:
onedayill wrote:
Car Dealership X, the total profit from sales increased by 10 percent over the previous years, while the number of cars sold decreased by 10 percent over the previous years. Approximately what was the average percent increase in profit per car over the previous year.

18%
20%
22%
23%
25%

Let Profit be x and cars sold be y

Originally Profit/Car = x/y

Now it is 1.1x / 0.9 y = 11/9 (x/y)

Increase in profit per car = ((2/9) (x/y)) / (x/y) * 100 = 22.22%

Thanks
I got my mistake.
_________________

GGG (Gym / GMAT / Girl) -- Be Serious

Its your duty to post OA afterwards; some one must be waiting for that...

Senior Manager
Joined: 20 Jul 2010
Posts: 263
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

23 Sep 2010, 20:52
I did this by putting numbers, I don't know if this approach is correct.

Previous Years - Number of cars be 10 , total profit be 100....

this year - Number of cars will be 9 , total profit will be 110..

Previous year ---- profit per car = 100/10 = 10
this year profit = 110/9 = 12.2
percent increase = 22%
_________________

If you like my post, consider giving me some KUDOS !!!!! Like you I need them

Retired Moderator
Joined: 02 Sep 2010
Posts: 803
Location: London
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

23 Sep 2010, 22:19
What is true in general will also be true in particular

Substituting numbers is perfectly fine
_________________
GMAT Club Legend
Joined: 09 Sep 2013
Posts: 16014
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

15 Dec 2015, 01:31
Hello from the GMAT Club BumpBot!

Thanks to another GMAT Club member, I have just discovered this valuable topic, yet it had no discussion for over a year. I am now bumping it up - doing my job. I think you may find it valuable (esp those replies with Kudos).

Want to see all other topics I dig out? Follow me (click follow button on profile). You will receive a summary of all topics I bump in your profile area as well as via email.
_________________
Director
Joined: 07 Dec 2014
Posts: 713
Re: Car Dealership X, the total profit from sales increased by 10 percent [#permalink]

### Show Tags

10 Apr 2016, 15:40
(1.1p/.9c)/(p/c)=1.22=22% increase
Re: Car Dealership X, the total profit from sales increased by 10 percent   [#permalink] 10 Apr 2016, 15:40
Similar topics Replies Last post
Similar
Topics:
1 The revenue from sales items in 1996 increases by x percent, ... 5 02 Oct 2016, 05:09
1 A house seller can deduct 5 percent from the total sales 5 25 Nov 2012, 20:17
1 A car dealership sold two cars: the first car at a 10% 1 17 Mar 2012, 04:33
9 A car dealership sold two cars: the first car at a 10% 16 06 Oct 2011, 08:10
23 Last September the sales total for a certain car dealership was 400% 9 30 Nov 2016, 14:16
Display posts from previous: Sort by

# Car Dealership X, the total profit from sales increased by 10 percent

 new topic post reply Question banks Downloads My Bookmarks Reviews Important topics

 Powered by phpBB © phpBB Group and phpBB SEO Kindly note that the GMAT® test is a registered trademark of the Graduate Management Admission Council®, and this site has neither been reviewed nor endorsed by GMAC®.


text:
# Eikonal Blog

## 2010.09.15

### Feynman derivation of Maxwell equations

• “On Feynman’s Approach to the Foundations of Gauge Theory” by M. C. Land, N. M. Shnerb and L. P. Horwitz; 2 Aug 1993. – 36 p. – http://cdsweb.cern.ch/record/568394/:
Abstract: In 1948, Feynman showed Dyson how the Lorentz force and Maxwell equations could be derived from commutation relations coordinates and velocities. Several authors noted that the derived equations are not Lorentz covariant and so are not the standard Maxwell theory. In particular, Hojman and Shepley proved that the existence of commutation relations is a strong assumption, sufficient to determine the corresponding action, which for Feynman’s derivation is of Newtonian form. Tanimura generalized Feynman’s derivation to a Lorentz covariant form, however, this derivation does not lead to the standard Maxwell theory either. Tanimura’s force equation depends on a fifth ({\it scalar}) electromagnetic potential, and the invariant evolution parameter cannot be consistently identified with the proper time of the particle motion. Moreover, the derivation cannot be made reparameterization invariant; the scalar potential causes violations of the mass-shell constraint which this invariance should guarantee. In this paper, we examine Tanimura’s derivation in the framework of the proper time method in relativistic mechanics, and use the technique of Hojman and Shepley to study the unconstrained commutation relations. We show that Tanimura’s result then corresponds to the five-dimensional electromagnetic theory previously derived from a Stueckelberg-type quantum theory in which one gauges the invariant parameter in the proper time method. This theory provides the final step in Feynman’s program of deriving the Maxwell theory from commutation relations; the Maxwell theory emerges as the “correlation limit” of a more general gauge theory, in which it is properly contained.

## 2010.02.05

### Equations of electro+magneto-dynamics in 3+1 dimensions

Filed under: physics — Tags: , — sandokan65 @ 22:49

Maxwell equations:

• $\nabla \vec{D} = \rho_e$,
• $\nabla \vec{B} = \rho_m$,
• $\partial_t \vec{B} + \nabla \times \vec{E} = - \vec{j}_m$,
• $\partial_t \vec{D} - \nabla \times \vec{H} = - \vec{j}_e$.

Integral form of Maxwell equations:

• $\oint_S \vec{D}\cdot d\vec{S} = Q_e :\equiv \int_V \rho_e dV$,
• $\oint_S \vec{B}\cdot d\vec{S} = Q_m :\equiv \int_V \rho_m dV$,
• $\oint_l \vec{E} \cdot d\vec{l} = - \int_S ( \vec{j}_m + \partial_t \vec{B})$
• $\oint_l \vec{H} \cdot d\vec{l} = + \int_S ( \vec{j}_e + \partial_t \vec{D})$

Moments of the field (aka Energy-Momentum(-Pressure) complex of the electromagnetic field):

• Flow of energy (Poynting vector): $\vec{\Gamma}:\equiv\vec{E}\times\vec{H}$,
• Linear momentum of field: $\vec{G}:\equiv \vec{B} \times\vec{D}$,
• Stress-pressure tensor: ${\cal T}:\equiv \frac12 \vec{E}\otimes\vec{D} + \frac12 \vec{H}\otimes\vec{B} - \frac14 {\bf 1} \frac12 (\vec{E}\cdot\vec{D} + \vec{H}\cdot\vec{B})$.

Material equations: in the case of the weak fields
the polarizable and magnetizable environment yields the linear responses:

• $\vec{D} = \epsilon \vec{E} + \vec{P} = \hat{\epsilon} \vec{E}$,
• $\vec{B} = \mu \vec{H} + \mu \vec{M} = \hat{\mu} \vec{H}$,

where the linear coefficients are the permittivity ($\epsilon$) and the the permeability ($\mu$), and the polarization and the magnetization of the substance are measured by vector fields $\vec{P}$ and $\vec{M}$.

When $\rho_m\ne 0$ or $\vec{\jmath}_m\ne 0$ one can not use the signgle complex of potentials $(\Phi, \vec{A})$ (as used in ordinary electrodynamics), but should add one more complex of potentials $(\Psi, \vec{F})$. Then:

• $\vec{E} = - \partial_t \vec{A} -\nabla \Phi -c^2 \nabla \times \vec{F}$,
• $\vec{D} = - \partial_t \vec{L} -\nabla \Psi - \nabla \times \vec{U}$,
• $\vec{H} = - \partial_t \vec{U} -\nabla \phi + c^2 \nabla \times \vec{L}$,
• $\vec{B} = - \partial_t \vec{F} -\nabla \psi + \nabla \times \vec{A}$.

The source constrainst:

• $\partial_t \rho_e + \nabla \cdot \vec{\jmath}_e = 0$,
• $\partial_t \rho_m + \nabla \cdot \vec{\jmath}_m = 0$.


text:
• Create Account

# Ellipsoid Collision Detection and Response - Incorrect Intersection Point

Old topic!

Guest, the last post of this topic is over 60 days old and at this point you may not reply in this topic. If you wish to continue this conversation start a new topic.

No replies to this topic

### #1BlueSpud  Members   -  Reputation: 543

Like
0Likes
Like

Posted 05 May 2014 - 05:56 PM

I've implemented the collision detection and response outlined in this paper:http://www.peroxide.dk/papers/collision/collision.pdf The collision detection of the whole thing works fine, but some quirky things happen when I try to do the collision response. Here is my current response code that (for the most part) works. - A side note: This code is just included in case the problem does not lie in the spots I have speculated it does. This code does not have to be read unless there is no obvious mistake in the code I have put further down below and have determined to be the culprit.

if (pysObjs.size() > 0)
{
//here we're going to test all the collision objects and apply the collision response to the nearest collision
glm::dvec3 velocity = playerESpace.vecInElipsoidSpace(glm::dvec3(-camx,-camy,-camz)-glm::dvec3(-lcamx,-lcamy,-lcamz));
glm::dvec3 basePos = playerESpace.vecInElipsoidSpace(glm::dvec3(-lcamx,-lcamy,-lcamz));
glm::dvec3 destination = playerESpace.vecInElipsoidSpace(glm::dvec3(-camx,-camy,-camz));
int object = 0;
double distance = 0.0;
bool initial = false;
std::vector<collisionResultsPacket>results;
for (int i = 0; i < pysObjs.size(); i++)
{
//we need to test each mesh and see which is the closest triangle from all the meshes
pysObjs[i].testCollisionMeshWithElipsoid(0);

//we just got the closest result for that mesh
//now we test if that is closer than the one we already have
if (initial == false || pysObjs[i].closestResults.distance < distance )
{
if (pysObjs[i].closestResults.packet.collided == true)
{
initial = true;
distance = pysObjs[i].closestResults.distance;
object = i;
}
}
}

//now we have the colsest collision point and we can calculate the collision response properly
if (pysObjs[object].closestResults.packet.collided == true )
{
basePos = basePos + velocity*(pysObjs[object].closestResults.packet.time0*0.3f);

glm::dvec3 newIntersectionPoint = pysObjs[object].closestResults.packet.intersectionPoint;

if (distance >= veryCloseDistance)
{
//velocity = glm::normalize(velocity)*(distance-veryCloseDistance);
//basePos = basePos+velocity;

velocity = glm::normalize(velocity);
newIntersectionPoint = pysObjs[object].closestResults.packet.intersectionPoint - (veryCloseDistance*velocity);

}

glm::dvec3 slidingNormal = glm::normalize(pysObjs[object].closestResults.packet.normal);

double distancetoCollision = glm::dot(destination, slidingNormal) - ((slidingNormal.x*basePos.x)+(slidingNormal.y*basePos.y)+(slidingNormal.z*basePos.z));
destination = destination - distancetoCollision*slidingNormal;

destination = -destination*glm::inverse(playerESpace.CBM);
basePos = -basePos*glm::inverse(playerESpace.CBM);

camx = destination.x;
camy = destination.y;
camz = destination.z;

lcamx = basePos.x;
lcamy = basePos.y;
lcamz = basePos.z;

}
}

You'll notice that the sliding plane's origin is not the intersection point.This is because the point is being calculated incorrectly as a bit of visual debugging has show. When the final position is set to the intersection point, the camera moves rather far away from the triangle. The code above works rather well, I just need to figure out why. I'v triple checked the code and I can't find any discrepancies between mine and the original, other than the pieces I've changed. Because none of the data in the response function is actually generated here, here is the code for the physicsObjects that actually check the collision.

//
// physicObject.cpp
// Injection
//
// Created by Logan on 8/7/13.
//

#include "physicObject.h"
{
model = ModelRegistry.accessModel(file);
collisionEnabled = true;
//we have to get the bounding box

}
bool physicsObject::onSameSide(glm::dvec3 point, glm::dvec3 p1, glm::dvec3 p2, glm::dvec3 p3)
{
glm::dvec3 cross1 = glm::cross(p3-p2, point-p2);
glm::dvec3 cross2 = glm::cross(p3-p2, p1-p2);
if (glm::dot(cross1, cross2) >= 0.0)
return true;
else return false;

}
collisionResponsePacket physicsObject::checkTriangle(unsigned int i)
{
collisionResponsePacket packetToReturn;
glm::dvec3 velocity = playerESpace.vecInElipsoidSpace(glm::dvec3(-camx,-camy,-camz)-glm::dvec3(-lcamx,-lcamy,-lcamz));
glm::dvec3 basePos = playerESpace.vecInElipsoidSpace(glm::dvec3(-lcamx,-lcamy,-lcamz));

//for now, we're just going to check the one triangle.
//.obj stores vertexes in counter-clockwise order, so we're going to swap the positions.
glm::dvec3 p1 = playerESpace.vecInElipsoidSpace((glm::dvec3(ModelRegistry.models[model].m.obj[i].x3,ModelRegistry.models[model].m.obj[i].y3,ModelRegistry.models[model].m.obj[i].z3)));
glm::dvec3 p2 = playerESpace.vecInElipsoidSpace((glm::dvec3(ModelRegistry.models[model].m.obj[i].x2,ModelRegistry.models[model].m.obj[i].y2,ModelRegistry.models[model].m.obj[i].z2)));
glm::dvec3 p3 = playerESpace.vecInElipsoidSpace((glm::dvec3(ModelRegistry.models[model].m.obj[i].x1,ModelRegistry.models[model].m.obj[i].y1,ModelRegistry.models[model].m.obj[i].z1)));

//get the normal
glm::dvec3 normal = glm::normalize(glm::cross(p2-p1, p3-p1));

//get signed distance
double distance = glm::dot(normal, basePos) - ((normal.x*p1.x)+(normal.y*p1.y)+(normal.z*p1.z));
double nvd = glm::dot(normal, velocity);

//if its equal to zero, the elipse is traveling parallel to the triangle
if (nvd != 0.0)
{
//figure out what time the sphere would intersect with the collision plane, it doesn't nessesarily have to
double time0,time1;
time0 = (-1.0f-distance)/nvd;
time1 = (1.0f-distance)/nvd;

if (time0 > time1)
{
double temp = time1;
time1 = time0;
time0 = temp;
}
//check if there can actually be a collision
if (time0 > 1.0f || time1 < 0.0f)
return packetToReturn;
else
{
//fix thew times if there was a possible collision
if (time0 < 0.0) time0 = 0.0;
if (time1 < 0.0) time1 = 0.0;
if (time0 > 1.0) time0 = 1.0;
if (time1 > 1.0) time1 = 1.0;

//continue to look for a collision
glm::dvec3 planeIntersectionPoint = (basePos-normal) + time0*velocity;
//now we have where the elipsoide will collide with the trangle, now we check if it is an actual collision
if (onSameSide(planeIntersectionPoint, p1, p2, p3) && onSameSide(planeIntersectionPoint, p2, p1, p3) && onSameSide(planeIntersectionPoint, p3, p1, p2))
{
//there was a collision here so we give the packet the proper information
packetToReturn.intersectionPoint = planeIntersectionPoint;
packetToReturn.time0 = time0;
packetToReturn.time1 = time1;
packetToReturn.collided = true;
packetToReturn.normal = normal;
}
else return packetToReturn;
}

} else return packetToReturn;

return packetToReturn;
}

void physicsObject::testCollisionMeshWithElipsoid(double time)
{
//some of the vectors we require for the response
glm::dvec3 velocity = playerESpace.vecInElipsoidSpace(glm::dvec3(-camx,-camy,-camz)-glm::dvec3(-lcamx,-lcamy,-lcamz));
std::vector<collisionResponsePacket>results;

for (int i = 0; i < ModelRegistry.models[model].m.obj.size(); i++)
//we chech each triangle
results.push_back(checkTriangle(i));

double closest = 0.0;
int triangle = 0;
bool initial = false;
//now we find which one we are going to respond against, then do it
for (int i = 0; i < results.size(); i++)
{
if (results[i].collided == true)
{
//get the distance to the collision
double distanceToCollision = results[i].time0*glm::length(velocity);
if (initial == false || distanceToCollision < closest )
{
//this one was less than the previos one
initial = true;
closest = distanceToCollision;
triangle = i;
}

}
}
//we store the information for the closest collision here
closestResults.i = triangle;
closestResults.distance = closest;
closestResults.packet = results[triangle];
results.clear();

}

The culprit I suspect is either in one of the two following places. This is the code I think needs to be looked at by a second pair of eyes, not necessarily the entire code above.

glm::dvec3 planeIntersectionPoint = (basePos-normal) + time0*velocity;

And the other piece I suspect could cause problems:

double distanceToCollision = results[i].time0*glm::length(velocity);

These are really the only things that are really used in the final collision detection. I'm satisfied with the outcome of the method that I implemented. Its just one of those things that keeps me awake at night. "Why didn't the original code work? Why didn't what the paper said work correctly?" The problem could easily lie in my implementation of the algorithm and if I had the choice, I'd re-write the whole program with this routine in mind. Unfortunately I can't do that. Any thought or input into why the original solution to the collision and response problem did not work for me would be much appreciated and it may help people having perhaps the same problem as me in the future.

Thank you.

Old topic!

Guest, the last post of this topic is over 60 days old and at this point you may not reply in this topic. If you wish to continue this conversation start a new topic.

PARTNERS


text:
# different types in same container

## Recommended Posts

sfx81    131
Hi. I was wondering is its possible to have a single array/container which can contain different types. e.g my array[]= { 1,'c', "string", customObject }; Regards

##### Share on other sites
popsoftheyear    2194
Advice: When asking questions like this, mention the language you are referring to. Otherwise, you're likely to get responses like "Well yeah, of course you can", which might be true for javascript, when you were actually talking about C++, where it isn't so true (well not naturally. But you could use a variant). Or whatever. I'm not even sure if javascript can do that, but you get my point.

So... what language?

Cheers
-Scott

##### Share on other sites
Antheus    2409
Yes, it's possible. But why would you want that?

Consider this: foo[17] returns something. What are you going to do with it? How will you know it's customObject and not a bool?

##### Share on other sites
sfx81    131
thanks every one for quick replies.
we are talking C++.

"foo[17] returns something. What are you going to do with it? How will you know it's customObject and not a bool?"
ANS: we can have a header for each object (say an enum) which tells what object it is, isnt it ?

I think varient is something I am after.
Cheers

##### Share on other sites
Telastyn    3777
Quote:
 Original post by sfx81ANS: we can have a header for each object (say an enum) which tells what object it is, isnt it ?

Which defeats the entire purpose of using a statically typed language.

##### Share on other sites
Antheus    2409
Quote:
 Original post by sfx81ANS: we can have a header for each object (say an enum) which tells what object it is, isnt it ?

This would typically be solved via polymorphism. Base class defines operations that can be performed on a type.
struct Base {  virtual void print() = 0;};struct IntValue : Base {  virtual void print() {    std::cout << value;  }  int value;};std::vector< Base* > values;for (int i = 0; i < values.size(); i++) {  values[i]->print();}

Quote:
 I think varient is something I am after.

Variant accomplishes the same as above, but instead of letting types define their own implementation of operations, visitor needs to know about all types. So visitor looks something like this:
struct Visitor {  void print(int i);  void print(std::string s);  void print(foo f);  ....};


text:
A Machine Learning Framework for Julia

# Introduction

MLJ (Machine Learning in Julia) is a toolbox written in Julia providing a common interface and meta-algorithms for selecting, tuning, evaluating, composing and comparing over 150 machine learning models written in Julia and other languages. In particular MLJ wraps a large number of scikit-learn models.

MLJ is released under the MIT licensed and sponsored by the Alan Turing Institute.

## Lightning tour

For a more elementary introduction to MLJ usage see Getting Started.

The first code snippet below creates a new Julia environment MLJ_tour and installs just those packages needed for the tour. See Installation for more on creating a Julia environment for use with MLJ.

Julia installation instructions are here.

using Pkg
Pkg.activate("MLJ_tour", shared=true)
Pkg.add("EvoTrees")

Load a selection of features and labels from the Ames House Price dataset:

using MLJ
X, y = @load_reduced_ames;

Load and instantiate a gradient tree-boosting model:

booster = @load EvoTreeRegressor
booster.max_depth = 2
booster.nrounds=50

Combine the model with categorical feature encoding:

pipe = @pipeline ContinuousEncoder booster

Define a hyper-parameter range for optimization:

max_depth_range = range(pipe,
:(evo_tree_regressor.max_depth),
lower = 1,
upper = 10)

Wrap the pipeline model in an optimization strategy:

self_tuning_pipe = TunedModel(model=pipe,
tuning=RandomSearch(),
ranges = max_depth_range,
resampling=CV(nfolds=3, rng=456),
measure=l1,
n=50)

Bind the "self-tuning" pipeline model (just a container for hyper-parameters) to data in a machine (which will additionally store learned parameters):

mach = machine(self_tuning_pipe, X, y)

Evaluate the "self-tuning" pipeline model's performance (implies nested resampling):

julia> evaluate!(mach,
measures=[l1, l2],
resampling=CV(nfolds=6, rng=123),
acceleration=CPUProcesses(), verbosity=2)
┌───────────┬───────────────┬────────────────────────────────────────────────────────┐
│ _.measure │ _.measurement │ _.per_fold                                             │
├───────────┼───────────────┼────────────────────────────────────────────────────────┤
│ l1        │ 16700.0       │ [16100.0, 16400.0, 14500.0, 17000.0, 16400.0, 19500.0] │
│ l2        │ 6.43e8        │ [5.88e8, 6.81e8, 4.35e8, 6.35e8, 5.98e8, 9.18e8]       │
└───────────┴───────────────┴────────────────────────────────────────────────────────┘
_.per_observation = [[[29100.0, 9990.0, ..., 103.0], [12100.0, 1330.0, ..., 13200.0], [6490.0, 22000.0, ..., 13800.0], [9090.0, 9530.0, ..., 13900.0], [50800.0, 22700.0, ..., 1550.0], [32800.0, 4940.0, ..., 1110.0]], [[8.45e8, 9.98e7, ..., 10500.0], [1.46e8, 1.77e6, ..., 1.73e8], [4.22e7, 4.86e8, ..., 1.9e8], [8.26e7, 9.09e7, ..., 1.93e8], [2.58e9, 5.13e8, ..., 2.42e6], [1.07e9, 2.44e7, ..., 1.24e6]]]
_.fitted_params_per_fold = [ … ]
_.report_per_fold = [ … ]

Try out MLJ yourself in the following batteries-included Binder notebook. No installation required.

## Key goals

• Offer a consistent way to use, compose and tune machine learning models in Julia,

• Promote the improvement of the Julia ML/Stats ecosystem by making it easier to use models from a wide range of packages,

• Unlock performance gains by exploiting Julia's support for parallelism, automatic differentiation, GPU, optimisation etc.

## Key features

• Data agnostic, train models on any data supported by the Tables.jl interface,

• Extensive, state-of-the art, support for model composition (pipelines and learning networks) (see more below),

• Convenient syntax to tune and evaluate (composite) models.

• Consistent interface to handle probabilistic predictions.

• Extensible tuning interface, to support growing number of optimization strategies, and designed to play well with model composition.

## Model composability

The generic model composition API's provided by other toolboxes we have surveyed share one or more of the following shortcomings, which do not exist in MLJ:

• Composite models do not inherit all the behavior of ordinary models.

• Composition is limited to linear (non-branching) pipelines.

• Supervised components in a linear pipeline can only occur at the end of the pipeline.

• Only static (unlearned) target transformations/inverse transformations are supported.

• Hyper-parameters in homogeneous model ensembles cannot be coupled.

• Model stacking, with out-of-sample predictions for base learners, cannot be implemented (using the generic API alone).

• Hyper-parameters and/or learned parameters of component models are not easily inspected or manipulated (by tuning algorithms, for example)

• Composite models cannot implement multiple opertations, for example, both a predict and transform method (as in clustering models) or both a transform and inverse_transform method.

Some of these features are demonstrated in this notebook

For more information see the MLJ design paper

## Reporting problems

Users are encouraged to provide feedback on their experience using MLJ and to report issues. You can do so here or on the #mlj Julia slack channel.

For known issues that are not strictly MLJ bugs, see here

## Installation

Initially it is recommended that MLJ and associated packages be installed in a new environment to avoid package conflicts. You can do this with

julia> using Pkg; Pkg.activate("my_MLJ_env", shared=true)

Installing MLJ is also done with the package manager:

julia> Pkg.add("MLJ")

Optional: To test your installation, run

julia> Pkg.test("MLJ")

It is important to note that MLJ is essentially a big wrapper providing a unified access to model providing packages. For this reason, one generally needs to add further packages to your environment to make model-specific code available. For instance, if you want to use a Decision Tree Classifier, you need to have MLJDecisionTreeInterface.jl installed:

julia> Pkg.add("MLJDecisionTreeInterface");
julia> using MLJ;
julia> @load DecisionTreeClassifier

However, if you try to use @load without adding the required package to your environment, an error message will tell you what package needs adding.

For a list of models and their packages run

using MLJ
models()

or refer to List of Supported Models

It is recommended that you start with models marked as coming from mature packages such as DecisionTree.jl, ScikitLearn.jl or XGBoost.jl.

MLJ is supported by a number of satelite packages (MLJTuning, MLJModelInterface, etc) which the general user is not required to install directly. Developers can learn more about these here

## Learning Julia

If you have experience in programming in another language but are new to Julia, then we highly recommend Aaron Christinson's tutorial Dispatching Design Patterns which is nicely compressed in his half-hour video presentation.

However, one doesn't need to be able to program in Julia to start using MLJ.

## Learning to use MLJ

The present document, although littered with examples, is primarily intended as a complete reference. For a short introduction to basic MLJ functionality, read the Getting Started section of this manual. For extensive tutorials, we recommend the MLJ Tutorials website. Each tutorial can be downloaded as a notebook or Julia script to facilitate experimentation. Finally, you may like to checkout the JuliaCon2020 Workshop on MLJ (recorded here).

You can try also MLJ out in the following notebook on Binder, without installing Julia or MLJ.

Users are also welcome to join the #mlj Julia slack channel to ask questions and make suggestions.

## Citing MLJ

When presenting work that uses MLJ, please cite the MLJ design paper:

@article{Blaom2020,
doi = {10.21105/joss.02704},
url = {https://doi.org/10.21105/joss.02704},
year = {2020},
publisher = {The Open Journal},
volume = {5},
number = {55},
pages = {2704},
author = {Anthony D. Blaom and Franz Kiraly and Thibaut Lienart and Yiannis Simillides and Diego Arenas and Sebastian J. Vollmer},
title = {{MLJ}: A Julia package for composable machine learning},
journal = {Journal of Open Source Software}
}


text:
## Chemistry and Chemical Reactivity (9th Edition)

(a) $Na_2C_2O_4$ is formed by: $Na^+$ (which doesn't affect the pH), and $C_2{O_4}^{2-}$ (which is a base); When we add a base to a solution, its pH is going to increase, because it is getting more basic. (b) Ammonium chloride: $NH_4Cl$: - $N{H_4}^+$ : An acid. - $Cl^-$: Neutral compound. When we add an acid to a solution, its pH is going to decrease, because it is getting more acidic. (c) $NaCl$ is formed by: $Na^+$ (which doesn't affect the pH) and $Cl^-$ (which also doesn't); Therefore, the pH will not change.


text:
## Results (displaying matches 1-50 of at least 1000) Next

Label Polynomial Discriminant Galois group Class group
1.1.1.1 x $$1$$ Trivial (as 1T1) Trivial
2.0.3.1 x2 - x + 1 $$-\,3$$ $C_2$ (as 2T1) Trivial
2.0.4.1 x2 + 1 $$-\,2^{2}$$ $C_2$ (as 2T1) Trivial
2.2.5.1 x2 - x - 1 $$5$$ $C_2$ (as 2T1) Trivial
2.0.7.1 x2 - x + 2 $$-\,7$$ $C_2$ (as 2T1) Trivial
2.0.8.1 x2 + 2 $$-\,2^{3}$$ $C_2$ (as 2T1) Trivial
2.2.8.1 x2 - 2 $$2^{3}$$ $C_2$ (as 2T1) Trivial
2.0.11.1 x2 - x + 3 $$-\,11$$ $C_2$ (as 2T1) Trivial
2.2.12.1 x2 - 3 $$2^{2}\cdot 3$$ $C_2$ (as 2T1) Trivial
2.2.13.1 x2 - x - 3 $$13$$ $C_2$ (as 2T1) Trivial
2.2.17.1 x2 - x - 4 $$17$$ $C_2$ (as 2T1) Trivial
2.0.19.1 x2 - x + 5 $$-\,19$$ $C_2$ (as 2T1) Trivial
2.2.21.1 x2 - x - 5 $$3\cdot 7$$ $C_2$ (as 2T1) Trivial
2.2.24.1 x2 - 6 $$2^{3}\cdot 3$$ $C_2$ (as 2T1) Trivial
2.2.28.1 x2 - 7 $$2^{2}\cdot 7$$ $C_2$ (as 2T1) Trivial
2.2.29.1 x2 - x - 7 $$29$$ $C_2$ (as 2T1) Trivial
2.2.33.1 x2 - x - 8 $$3\cdot 11$$ $C_2$ (as 2T1) Trivial
2.2.37.1 x2 - x - 9 $$37$$ $C_2$ (as 2T1) Trivial
2.2.41.1 x2 - x - 10 $$41$$ $C_2$ (as 2T1) Trivial
2.0.43.1 x2 - x + 11 $$-\,43$$ $C_2$ (as 2T1) Trivial
2.2.44.1 x2 - 11 $$2^{2}\cdot 11$$ $C_2$ (as 2T1) Trivial
2.2.53.1 x2 - x - 13 $$53$$ $C_2$ (as 2T1) Trivial
2.2.56.1 x2 - 14 $$2^{3}\cdot 7$$ $C_2$ (as 2T1) Trivial
2.2.57.1 x2 - x - 14 $$3\cdot 19$$ $C_2$ (as 2T1) Trivial
2.2.61.1 x2 - x - 15 $$61$$ $C_2$ (as 2T1) Trivial
2.0.67.1 x2 - x + 17 $$-\,67$$ $C_2$ (as 2T1) Trivial
2.2.69.1 x2 - x - 17 $$3\cdot 23$$ $C_2$ (as 2T1) Trivial
2.2.73.1 x2 - x - 18 $$73$$ $C_2$ (as 2T1) Trivial
2.2.76.1 x2 - 19 $$2^{2}\cdot 19$$ $C_2$ (as 2T1) Trivial
2.2.77.1 x2 - x - 19 $$7\cdot 11$$ $C_2$ (as 2T1) Trivial
2.2.88.1 x2 - 22 $$2^{3}\cdot 11$$ $C_2$ (as 2T1) Trivial
2.2.89.1 x2 - x - 22 $$89$$ $C_2$ (as 2T1) Trivial
2.2.92.1 x2 - 23 $$2^{2}\cdot 23$$ $C_2$ (as 2T1) Trivial
2.2.93.1 x2 - x - 23 $$3\cdot 31$$ $C_2$ (as 2T1) Trivial
2.2.97.1 x2 - x - 24 $$97$$ $C_2$ (as 2T1) Trivial
2.2.101.1 x2 - x - 25 $$101$$ $C_2$ (as 2T1) Trivial
2.2.109.1 x2 - x - 27 $$109$$ $C_2$ (as 2T1) Trivial
2.2.113.1 x2 - x - 28 $$113$$ $C_2$ (as 2T1) Trivial
2.2.124.1 x2 - 31 $$2^{2}\cdot 31$$ $C_2$ (as 2T1) Trivial
2.2.129.1 x2 - x - 32 $$3\cdot 43$$ $C_2$ (as 2T1) Trivial
2.2.133.1 x2 - x - 33 $$7\cdot 19$$ $C_2$ (as 2T1) Trivial
2.2.137.1 x2 - x - 34 $$137$$ $C_2$ (as 2T1) Trivial
2.2.141.1 x2 - x - 35 $$3\cdot 47$$ $C_2$ (as 2T1) Trivial
2.2.149.1 x2 - x - 37 $$149$$ $C_2$ (as 2T1) Trivial
2.2.152.1 x2 - 38 $$2^{3}\cdot 19$$ $C_2$ (as 2T1) Trivial
2.2.157.1 x2 - x - 39 $$157$$ $C_2$ (as 2T1) Trivial
2.2.161.1 x2 - x - 40 $$7\cdot 23$$ $C_2$ (as 2T1) Trivial
2.0.163.1 x2 - x + 41 $$-\,163$$ $C_2$ (as 2T1) Trivial
2.2.172.1 x2 - 43 $$2^{2}\cdot 43$$ $C_2$ (as 2T1) Trivial
2.2.173.1 x2 - x - 43 $$173$$ $C_2$ (as 2T1) Trivial
Next


text:
# Create Title Page with Hepthesis produces error

I'm trying to create a title page with hepthesis but run into problems. Could you hint me in the direction where I skrewed up?

MWE1:

\documentclass[oneside,a4paper,]{hepthesis}

\begin{document}
\author{QWER QWER}
\title{ABC acb}
\date{1.1.2000}

\begin{titlepage}
BlaBla
\end{titlepage}

\end{document}

procudes: ! LaTeX Error: There's no line here to end.

MWE2:

\documentclass[oneside,a4paper,]{hepthesis}

\begin{document}
\author{QWER QWER}
\title{ABC acb}
\date{1.1.2000}

\begin{titlepage}
\begin{minipage}{0.2\textwidth}
asdf
\end{minipage}
\end{titlepage}

\end{document}

produces the same with a lot of additional errors which migth go away once I takle the problem. Omitting the lines of author/title/date didn't help.

-

In hepthesis, there is no titlepage environment, but a \titlepage command, see section 6.5 in the class documentation.

The optional command is not really optional, because its absence will still cause the There's no line here to end, telling me that the class has not really been fully tested. Also, an input such as

\documentclass[oneside,a4paper,]{hepthesis}

\begin{document}
\author{QWER QWER}
\title{ABC acb}
\date{1.1.2000}

\titlepage[Some affiliation]{%
A dissertation submitted to the University of Somewhere\\
for the degree of Doctor of Philosophy
}

\end{document}

will produce the following typographically disputable result:

A sounder definition of \titlepage would be

\documentclass[oneside,a4paper,]{hepthesis}

\makeatletter
\renewcommand{\titlepage}[2][]{%
\par
\begingroup
\thispagestyle{empty}%
\ifx\@sftitles\@empty\else\sffamily\fi%
\centering
\vspace*{\frontmattertitleskip}%
\begin{doublespace}%
\Huge\textbf{\thetitle}\\
\end{doublespace}%
\vspace*{3cm}%
{\Large\theauthor\\[1ex]}
{#1\par}
\vspace*{8cm}%
{#2\\}
\endgroup
}
\makeatother

\begin{document}
\author{QWER QWER}
\title{ABC acb}
\date{1.1.2000}

\titlepage[Some affiliation]{%
A dissertation submitted to the University of Somewhere\\
for the degree of Doctor of Philosophy
}

\end{document}

This will work also without the optional argument.

-
Okey Thank you for your answer, I'll use the renewcommand as my titlepage environment now. Seems a bit hacky but this might work. –  magu_ Jul 18 '14 at 14:32


text:
The sum of length, breadth and depth of a cuboid is 19 cm and the length of its diagonal is 11 cm.
Question:

The sum of length, breadth and depth of a cuboid is 19 cm and the length of its diagonal is 11 cm. Find the surface area of the cuboid.

Solution:

Let the length, breadth and height (or depth) of the cuboid be l cm, b cm and h cm, respectively.

∴ l + b + h = 19           …..(1)

Also,

Length of the diagonal = 11 cm

$\Rightarrow \sqrt{l^{2}+b^{2}+h^{2}}=11$

$\Rightarrow l^{2}+b^{2}+h^{2}=121$           ……(2)

Squaring (1), we get

(l + b + h)2 = 192

⇒ l2 + b2 + h+ 2(lb + bh + hl) = 361

⇒ 121 + 2(lb + bh + hl) = 361                      [Using (2)]

⇒ 2(lb + bh + hl) = 361 − 121 = 240 cm2

Thus, the surface area of the cuboid is 240 cm2.