text:
A very wide range of physical processes lead to wave motion, where signals are propagated through a medium in space and time, normally with little or no permanent movement of the medium itself. The shape of the signals may undergo changes as they travel through matter, but usually not so much that the signals cannot be recognized at some later point in space and time. Many types of wave motion can be described by the equation $$u_{tt}=\nabla\cdot (c^2\nabla u) + f$$, which we will solve in the forthcoming text by finite difference methods.

# Simulation of waves on a string¶

We begin our study of wave equations by simulating one-dimensional waves on a string, say on a guitar or violin. Let the string in the undeformed state coincide with the interval $$[0,L]$$ on the $$x$$ axis, and let $$u(x,t)$$ be the displacement at time $$t$$ in the $$y$$ direction of a point initially at $$x$$. The displacement function $$u$$ is governed by the mathematical model

$$$\frac{\partial^2 u}{\partial t^2} = c^2 \frac{\partial^2 u}{\partial x^2}, \quad x\in (0,L),\ t\in (0,T] \label{wave:pde1} \tag{1}$$$
$$$u(x,0) = I(x), \quad x\in [0,L] \label{wave:pde1:ic:u} \tag{2}$$$
$$$\frac{\partial}{\partial t}u(x,0) = 0, \quad x\in [0,L] \label{wave:pde1:ic:ut} \tag{3}$$$
$$$u(0,t) = 0, \quad t\in (0,T] \label{wave:pde1:bc:0} \tag{4}$$$
$$$u(L,t) = 0, \quad t\in (0,T] \label{wave:pde1:bc:L} \tag{5}$$$

The constant $$c$$ and the function $$I(x)$$ must be prescribed.

Equation (1) is known as the one-dimensional wave equation. Since this PDE contains a second-order derivative in time, we need two initial conditions. The condition (2) specifies the initial shape of the string, $$I(x)$$, and (3) expresses that the initial velocity of the string is zero. In addition, PDEs need boundary conditions, given here as (4) and (5). These two conditions specify that the string is fixed at the ends, i.e., that the displacement $$u$$ is zero.

The solution $$u(x,t)$$ varies in space and time and describes waves that move with velocity $$c$$ to the left and right.

Sometimes we will use a more compact notation for the partial derivatives to save space:

$$$u_t = \frac{\partial u}{\partial t}, \quad u_{tt} = \frac{\partial^2 u}{\partial t^2}, \label{_auto1} \tag{6}$$$

and similar expressions for derivatives with respect to other variables. Then the wave equation can be written compactly as $$u_{tt} = c^2u_{xx}$$.

The PDE problem (1)-(5) will now be discretized in space and time by a finite difference method.

## Discretizing the domain¶

The temporal domain $$[0,T]$$ is represented by a finite number of mesh points

$$$0 = t_0 < t_1 < t_2 < \cdots < t_{N_t-1} < t_{N_t} = T \label{_auto2} \tag{7}$$$

Similarly, the spatial domain $$[0,L]$$ is replaced by a set of mesh points

$$$0 = x_0 < x_1 < x_2 < \cdots < x_{N_x-1} < x_{N_x} = L \label{_auto3} \tag{8}$$$

One may view the mesh as two-dimensional in the $$x,t$$ plane, consisting of points $$(x_i, t_n)$$, with $$i=0,\ldots,N_x$$ and $$n=0,\ldots,N_t$$.

### Uniform meshes¶

For uniformly distributed mesh points we can introduce the constant mesh spacings $$\Delta t$$ and $$\Delta x$$. We have that

$$$x_i = i\Delta x,\ i=0,\ldots,N_x,\quad t_n = n\Delta t,\ n=0,\ldots,N_t \label{_auto4} \tag{9}$$$

We also have that $$\Delta x = x_i-x_{i-1}$$, $$i=1,\ldots,N_x$$, and $$\Delta t = t_n - t_{n-1}$$, $$n=1,\ldots,N_t$$. Figure displays a mesh in the $$x,t$$ plane with $$N_t=5$$, $$N_x=5$$, and constant mesh spacings.

## The discrete solution¶

The solution $$u(x,t)$$ is sought at the mesh points. We introduce the mesh function $$u_i^n$$, which approximates the exact solution at the mesh point $$(x_i,t_n)$$ for $$i=0,\ldots,N_x$$ and $$n=0,\ldots,N_t$$. Using the finite difference method, we shall develop algebraic equations for computing the mesh function.

## Fulfilling the equation at the mesh points¶

In the finite difference method, we relax the condition that (1) holds at all points in the space-time domain $$(0,L)\times (0,T]$$ to the requirement that the PDE is fulfilled at the interior mesh points only:

$$$\frac{\partial^2}{\partial t^2} u(x_i, t_n) = c^2\frac{\partial^2}{\partial x^2} u(x_i, t_n), \label{wave:pde1:step2} \tag{10}$$$

for $$i=1,\ldots,N_x-1$$ and $$n=1,\ldots,N_t-1$$. For $$n=0$$ we have the initial conditions $$u=I(x)$$ and $$u_t=0$$, and at the boundaries $$i=0,N_x$$ we have the boundary condition $$u=0$$.

## Replacing derivatives by finite differences¶

The second-order derivatives can be replaced by central differences. The most widely used difference approximation of the second-order derivative is

$\frac{\partial^2}{\partial t^2}u(x_i,t_n)\approx \frac{u_i^{n+1} - 2u_i^n + u^{n-1}_i}{\Delta t^2}\$

mathcal{I}_t is convenient to introduce the finite difference operator notation

$[D_tD_t u]^n_i = \frac{u_i^{n+1} - 2u_i^n + u^{n-1}_i}{\Delta t^2}$

A similar approximation of the second-order derivative in the $$x$$ direction reads

$\frac{\partial^2}{\partial x^2}u(x_i,t_n)\approx \frac{u_{i+1}^{n} - 2u_i^n + u^{n}_{i-1}}{\Delta x^2} = [D_xD_x u]^n_i$

### Algebraic version of the PDE¶

We can now replace the derivatives in (10) and get

$$$\frac{u_i^{n+1} - 2u_i^n + u^{n-1}_i}{\Delta t^2} = c^2\frac{u_{i+1}^{n} - 2u_i^n + u^{n}_{i-1}}{\Delta x^2}, \label{wave:pde1:step3b} \tag{11}$$$

or written more compactly using the operator notation:

$$$[D_tD_t u = c^2 D_xD_x]^{n}_i \label{wave:pde1:step3a} \tag{12}$$$

### Interpretation of the equation as a stencil¶

A characteristic feature of (11) is that it involves $$u$$ values from neighboring points only: $$u_i^{n+1}$$, $$u^n_{i\pm 1}$$, $$u^n_i$$, and $$u^{n-1}_i$$. The circles in Figure illustrate such neighboring mesh points that contribute to an algebraic equation. In this particular case, we have sampled the PDE at the point $$(2,2)$$ and constructed (11), which then involves a coupling of $$u_1^2$$, $$u_2^3$$, $$u_2^2$$, $$u_2^1$$, and $$u_3^2$$. The term stencil is often used about the algebraic equation at a mesh point, and the geometry of a typical stencil is illustrated in Figure. One also often refers to the algebraic equations as discrete equations, (finite) difference equations or a finite difference scheme.

Mesh in space and time. The circles show points connected in a finite difference equation.

### Algebraic version of the initial conditions¶

We also need to replace the derivative in the initial condition (3) by a finite difference approximation. A centered difference of the type

$\frac{\partial}{\partial t} u(x_i,t_0)\approx \frac{u^1_i - u^{-1}_i}{2\Delta t} = [D_{2t} u]^0_i,$

seems appropriate. Writing out this equation and ordering the terms give

$$$u^{-1}_i=u^{1}_i,\quad i=0,\ldots,N_x \label{wave:pde1:step3c} \tag{13}$$$

The other initial condition can be computed by

$u_i^0 = I(x_i),\quad i=0,\ldots,N_x$

## Formulating a recursive algorithm¶

We assume that $$u^n_i$$ and $$u^{n-1}_i$$ are available for $$i=0,\ldots,N_x$$. The only unknown quantity in (11) is therefore $$u^{n+1}_i$$, which we now can solve for:

$$$u^{n+1}_i = -u^{n-1}_i + 2u^n_i + C^2 \left(u^{n}_{i+1}-2u^{n}_{i} + u^{n}_{i-1}\right) \label{wave:pde1:step4} \tag{14}$$$

We have here introduced the parameter

$$$C = c\frac{\Delta t}{\Delta x}, \label{_auto5} \tag{15}$$$

known as the Courant number.

$$C$$ is the key parameter in the discrete wave equation.

We see that the discrete version of the PDE features only one parameter, $$C$$, which is therefore the key parameter, together with $$N_x$$, that governs the quality of the numerical solution (see the section Analysis of the difference equations for details). Both the primary physical parameter $$c$$ and the numerical parameters $$\Delta x$$ and $$\Delta t$$ are lumped together in $$C$$. Note that $$C$$ is a dimensionless parameter.

Given that $$u^{n-1}_i$$ and $$u^n_i$$ are known for $$i=0,\ldots,N_x$$, we find new values at the next time level by applying the formula (14) for $$i=1,\ldots,N_x-1$$. Figure illustrates the points that are used to compute $$u^3_2$$. For the boundary points, $$i=0$$ and $$i=N_x$$, we apply the boundary conditions $$u_i^{n+1}=0$$.

Even though sound reasoning leads up to (14), there is still a minor challenge with it that needs to be resolved. Think of the very first computational step to be made. The scheme (14) is supposed to start at $$n=1$$, which means that we compute $$u^2$$ from $$u^1$$ and $$u^0$$. Unfortunately, we do not know the value of $$u^1$$, so how to proceed? A standard procedure in such cases is to apply (14) also for $$n=0$$. This immediately seems strange, since it involves $$u^{-1}_i$$, which is an undefined quantity outside the time mesh (and the time domain). However, we can use the initial condition (13) in combination with (14) when $$n=0$$ to eliminate $$u^{-1}_i$$ and arrive at a special formula for $$u_i^1$$:

$$$u_i^1 = u^0_i - \frac{1}{2} C^2\left(u^{0}_{i+1}-2u^{0}_{i} + u^{0}_{i-1}\right) \label{wave:pde1:step4:1} \tag{16}$$$

Figure illustrates how (16) connects four instead of five points: $$u^1_2$$, $$u_1^0$$, $$u_2^0$$, and $$u_3^0$$.

Modified stencil for the first time step.

We can now summarize the computational algorithm:

1. Compute $$u^0_i=I(x_i)$$ for $$i=0,\ldots,N_x$$

2. Compute $$u^1_i$$ by (16) for $$i=1,2,\ldots,N_x-1$$ and set $$u_i^1=0$$ for the boundary points given by $$i=0$$ and $$i=N_x$$,

3. For each time level $$n=1,2,\ldots,N_t-1$$

a. apply (14) to find $$u^{n+1}_i$$ for $$i=1,\ldots,N_x-1$$

b. set $$u^{n+1}_i=0$$ for the boundary points having $$i=0$$, $$i=N_x$$.

The algorithm essentially consists of moving a finite difference stencil through all the mesh points, which can be seen as an animation in a web page or a movie file.

## Sketch of an implementation¶

We start by defining some constants that will be used throughout our Devito code.

import numpy as np

# Given mesh points as arrays x and t (x[i], t[n]),
# constant c and function I for initial condition
x = np.linspace(0, 2, 101)
t = np.linspace(0, 2, 101)
c = 1
I = lambda x: np.sin(x)

dx = x[1] - x[0]
dt = t[1] - t[0]
C = c*dt/dx            # Courant number
Nx = len(x)-1
Nt = len(t)-1
C2 = C**2              # Help variable in the scheme
L = 2.

Next, we define our 1D computational grid and create a function u as a symbolic devito.TimeFunction. We need to specify the space_order as 2 since our wave equation involves second-order derivatives with respect to $$x$$. Similarly, we specify the time_order as 2, as our equation involves second-order derivatives with respect to $$t$$. Setting these parameters allows us to use u.dx2 and u.dt2.

from devito import Grid, TimeFunction

# Initialise u for space and time order 2, using initialisation function I
grid = Grid(shape=(Nx+1), extent=(L))
u = TimeFunction(name='u', grid=grid, time_order=2, space_order=2)
u.data[:,:] = I(x[:])

Now that we have initialised u, we can solve our wave equation for the unknown quantity $$u^{n+1}_i$$ using forward and backward differences in space and time.

from devito import Constant, Eq, solve

# Set up wave equation and solve for forward stencil point in time
pde = (1/c**2)*u.dt2-u.dx2
stencil = Eq(u.forward, solve(pde, u.forward))

print("LHS: %s" % stencil.lhs)
print("RHS: %s" % stencil.rhs)

LHS: u(t + dt, x)
RHS: 1.0*dt**2*(-2.0*u(t, x)/h_x**2 + u(t, x - h_x)/h_x**2 + u(t, x + h_x)/h_x**2 + 2.0*u(t, x)/dt**2 - 1.0*u(t - dt, x)/dt**2)

Great! From these print statements, we can see that Devito has taken the wave equation in (1) and solved it for $$u^{n+1}_i$$, giving us equation (14). Note that dx is denoted as h_x, while u(t, x), u(t, x - h_x) and u(t, x + h_x) denote the equivalent of $$u^{n}_{i}$$, $$u^{n}_{i-1}$$ and $$u^{n}_{i+1}$$ respectively.

We also need to create a separate stencil for the first timestep, where we substitute $$u^{1}_i$$ for $$u^{-1}_i$$, as given in (13).

stencil_init = stencil.subs(u.backward, u.forward)

Now we can create expressions for our boundary conditions and build the operator. The results are plotted below.

#NBVAL_IGNORE_OUTPUT
from devito import Operator

t_s = grid.stepping_dim

# Boundary conditions
bc = [Eq(u[t_s+1, 0], 0)]
bc += [Eq(u[t_s+1, Nx], 0)]

# Defining one Operator for initial timestep and one for the rest
op_init = Operator([stencil_init]+bc)
op = Operator([stencil]+bc)

op_init.apply(time_M=1, dt=dt)
op.apply(time_m=1,time_M=Nt, dt=dt)

Data type float64 of runtime value dt does not match the Constant data type <class 'numpy.float32'>

Operator Kernel run in 0.01 s

Data type float64 of runtime value dt does not match the Constant data type <class 'numpy.float32'>

Operator Kernel run in 0.01 s

PerformanceSummary([(PerfKey(name='section0', rank=None),
PerfEntry(time=2.8000000000000003e-05, gflopss=0.0, gpointss=0.0, oi=0.0, ops=0, itershapes=[]))])

We can plot our results using matplotlib:

import matplotlib.pyplot as plt

plt.plot(x, u.data[-1])
plt.xlabel('x')
plt.ylabel('u')
plt.show()

# Verification¶

Before implementing the algorithm, it is convenient to add a source term to the PDE (1), since that gives us more freedom in finding test problems for verification. Physically, a source term acts as a generator for waves in the interior of the domain.

## A slightly generalized model problem¶

We now address the following extended initial-boundary value problem for one-dimensional wave phenomena:

$$$u_{tt} = c^2 u_{xx} + f(x,t), \quad x\in (0,L),\ t\in (0,T] \label{wave:pde2} \tag{17}$$$
$$$u(x,0) = I(x), \quad x\in [0,L] \label{wave:pde2:ic:u} \tag{18}$$$
$$$u_t(x,0) = V(x), \quad x\in [0,L] \label{wave:pde2:ic:ut} \tag{19}$$$
$$$u(0,t) = 0, \quad t>0 \label{wave:pde2:bc:0} \tag{20}$$$
$$$u(L,t) = 0, \quad t>0 \label{wave:pde2:bc:L} \tag{21}$$$

Sampling the PDE at $$(x_i,t_n)$$ and using the same finite difference approximations as above, yields

$$$[D_tD_t u = c^2 D_xD_x u + f]^{n}_i \label{wave:pde2:fdop} \tag{22}$$$

Writing this out and solving for the unknown $$u^{n+1}_i$$ results in

$$$u^{n+1}_i = -u^{n-1}_i + 2u^n_i + C^2 (u^{n}_{i+1}-2u^{n}_{i} + u^{n}_{i-1}) + \Delta t^2 f^n_i \label{wave:pde2:step3b} \tag{23}$$$

The equation for the first time step must be rederived. The discretization of the initial condition $$u_t = V(x)$$ at $$t=0$$ becomes

$[D_{2t}u = V]^0_i\quad\Rightarrow\quad u^{-1}_i = u^{1}_i - 2\Delta t V_i,$

which, when inserted in (23) for $$n=0$$, gives the special formula

$$$u^{1}_i = u^0_i + \Delta t V_i + {\frac{1}{2}} C^2 \left(u^{0}_{i+1}-2u^{0}_{i} + u^{0}_{i-1}\right) + \frac{1}{2}\Delta t^2 f^0_i \label{wave:pde2:step3c} \tag{24}$$$

## Using an analytical solution of physical significance¶

Many wave problems feature sinusoidal oscillations in time and space. For example, the original PDE problem (1)-(5) allows an exact solution

$$$u_e(x,t) = A\sin\left(\frac{\pi}{L}x\right) \cos\left(\frac{\pi}{L}ct\right) \label{wave:pde2:test:ue} \tag{25}$$$

This $$u_e$$ fulfills the PDE with $$f=0$$, boundary conditions $$u_e(0,t)=u_e(L,t)=0$$, as well as initial conditions $$I(x)=A\sin\left(\frac{\pi}{L}x\right)$$ and $$V=0$$.

How to use exact solutions for verification.

It is common to use such exact solutions of physical interest to verify implementations. However, the numerical solution $$u^n_i$$ will only be an approximation to $$u_e(x_i,t_n)$$. We have no knowledge of the precise size of the error in this approximation, and therefore we can never know if discrepancies between $$u^n_i$$ and $$u_e(x_i,t_n)$$ are caused by mathematical approximations or programming errors. In particular, if plots of the computed solution $$u^n_i$$ and the exact one (25) look similar, many are tempted to claim that the implementation works. However, even if color plots look nice and the accuracy is “deemed good”, there can still be serious programming errors present!

The only way to use exact physical solutions like (25) for serious and thorough verification is to run a series of simulations on finer and finer meshes, measure the integrated error in each mesh, and from this information estimate the empirical convergence rate of the method.

An introduction to the computing of convergence rates is given in Section 3.1.6 in [Langtangen_decay]. There is also a detailed example on computing convergence rates in the verification section of the Vibration ODEs chapter.

In the present problem, one expects the method to have a convergence rate of 2 (see the section Analysis of the difference equations), so if the computed rates are close to 2 on a sufficiently fine mesh, we have good evidence that the implementation is free of programming mistakes.

## Manufactured solution and estimation of convergence rates¶

### Specifying the solution and computing corresponding data¶

One problem with the exact solution (25) is that it requires a simplification ($${V}=0, f=0$$) of the implemented problem (17)-(21). An advantage of using a manufactured solution is that we can test all terms in the PDE problem. The idea of this approach is to set up some chosen solution and fit the source term, boundary conditions, and initial conditions to be compatible with the chosen solution. Given that our boundary conditions in the implementation are $$u(0,t)=u(L,t)=0$$, we must choose a solution that fulfills these conditions. One example is

$u_e(x,t) = x(L-x)\sin t$

Inserted in the PDE $$u_{tt}=c^2u_{xx}+f$$ we get

$-x(L-x)\sin t = -c^2 2\sin t + f\quad\Rightarrow f = (2c^2 - x(L-x))\sin t$

The initial conditions become

\begin{split} \begin{align*} u(x,0) =& I(x) = 0,\\ u_t(x,0) &= V(x) = x(L-x) \end{align*} \end{split}

### Defining a single discretization parameter¶

To verify the code, we compute the convergence rates in a series of simulations, letting each simulation use a finer mesh than the previous one. Such empirical estimation of convergence rates relies on an assumption that some measure $$E$$ of the numerical error is related to the discretization parameters through

$E = C_t\Delta t^r + C_x\Delta x^p,$

where $$C_t$$, $$C_x$$, $$r$$, and $$p$$ are constants. The constants $$r$$ and $$p$$ are known as the convergence rates in time and space, respectively. From the accuracy in the finite difference approximations, we expect $$r=p=2$$, since the error terms are of order $$\Delta t^2$$ and $$\Delta x^2$$. This is confirmed by truncation error analysis and other types of analysis.

By using an exact solution of the PDE problem, we will next compute the error measure $$E$$ on a sequence of refined meshes and see if the rates $$r=p=2$$ are obtained. We will not be concerned with estimating the constants $$C_t$$ and $$C_x$$, simply because we are not interested in their values.

mathcal{I}_t is advantageous to introduce a single discretization parameter $$h=\Delta t=\hat c \Delta x$$ for some constant $$\hat c$$. Since $$\Delta t$$ and $$\Delta x$$ are related through the Courant number, $$\Delta t = C\Delta x/c$$, we set $$h=\Delta t$$, and then $$\Delta x = hc/C$$. Now the expression for the error measure is greatly simplified:

$E = C_t\Delta t^r + C_x\Delta x^r = C_t h^r + C_x\left(\frac{c}{C}\right)^r h^r = Dh^r,\quad D = C_t+C_x\left(\frac{c}{C}\right)^r$

### Computing errors¶

We choose an initial discretization parameter $$h_0$$ and run experiments with decreasing $$h$$: $$h_i=2^{-i}h_0$$, $$i=1,2,\ldots,m$$. Halving $$h$$ in each experiment is not necessary, but it is a common choice. For each experiment we must record $$E$$ and $$h$$. Standard choices of error measure are the $$\ell^2$$ and $$\ell^\infty$$ norms of the error mesh function $$e^n_i$$:

$$$E = ||e^n_i||_{\ell^2} = \left( \Delta t\Delta x \sum_{n=0}^{N_t}\sum_{i=0}^{N_x} (e^n_i)^2\right)^{\frac{1}{2}},\quad e^n_i = u_e(x_i,t_n)-u^n_i, \label{wave:pde2:fd:MMS:E:l2} \tag{26}$$$
$$$E = ||e^n_i||_{\ell^\infty} = \max_{i,n} |e^n_i| \label{wave:pde2:fd:MMS:E:linf} \tag{27}$$$

In Python, one can compute $$\sum_{i}(e^{n}_i)^2$$ at each time step and accumulate the value in some sum variable, say e2_sum. At the final time step one can do sqrt(dt*dx*e2_sum). For the $$\ell^\infty$$ norm one must compare the maximum error at a time level (e.max()) with the global maximum over the time domain: e_max = max(e_max, e.max()).

An alternative error measure is to use a spatial norm at one time step only, e.g., the end time $$T$$ ($$n=N_t$$):

$$$E = ||e^n_i||_{\ell^2} = \left( \Delta x\sum_{i=0}^{N_x} (e^n_i)^2\right)^{\frac{1}{2}},\quad e^n_i = u_e(x_i,t_n)-u^n_i, \label{_auto6} \tag{28}$$$
$$$E = ||e^n_i||_{\ell^\infty} = \max_{0\leq i\leq N_x} |e^{n}_i| \label{_auto7} \tag{29}$$$

The important point is that the error measure ($$E$$) for the simulation is represented by a single number.

### Computing rates¶

Let $$E_i$$ be the error measure in experiment (mesh) number $$i$$ (not to be confused with the spatial index $$i$$) and let $$h_i$$ be the corresponding discretization parameter ($$h$$). With the error model $$E_i = Dh_i^r$$, we can estimate $$r$$ by comparing two consecutive experiments:

\begin{split} \begin{align*} E_{i+1}& =D h_{i+1}^{r},\\ E_{i}& =D h_{i}^{r} \end{align*} \end{split}

Dividing the two equations eliminates the (uninteresting) constant $$D$$. Thereafter, solving for $$r$$ yields

$r = \frac{\ln E_{i+1}/E_{i}}{\ln h_{i+1}/h_{i}}$

Since $$r$$ depends on $$i$$, i.e., which simulations we compare, we add an index to $$r$$: $$r_i$$, where $$i=0,\ldots,m-2$$, if we have $$m$$ experiments: $$(h_0,E_0),\ldots,(h_{m-1}, E_{m-1})$$.

In our present discretization of the wave equation we expect $$r=2$$, and hence the $$r_i$$ values should converge to 2 as $$i$$ increases.

## Constructing an exact solution of the discrete equations¶

With a manufactured or known analytical solution, as outlined above, we can estimate convergence rates and see if they have the correct asymptotic behavior. Experience shows that this is a quite good verification technique in that many common bugs will destroy the convergence rates. A significantly better test though, would be to check that the numerical solution is exactly what it should be. This will in general require exact knowledge of the numerical error, which we do not normally have (although we in the section Analysis of the difference equations establish such knowledge in simple cases). However, it is possible to look for solutions where we can show that the numerical error vanishes, i.e., the solution of the original continuous PDE problem is also a solution of the discrete equations. This property often arises if the exact solution of the PDE is a lower-order polynomial. (Truncation error analysis leads to error measures that involve derivatives of the exact solution. In the present problem, the truncation error involves 4th-order derivatives of $$u$$ in space and time. Choosing $$u$$ as a polynomial of degree three or less will therefore lead to vanishing error.)

We shall now illustrate the construction of an exact solution to both the PDE itself and the discrete equations. Our chosen manufactured solution is quadratic in space and linear in time. More specifically, we set

$$$u_e (x,t) = x(L-x)(1+{\frac{1}{2}}t), \label{wave:pde2:fd:verify:quadratic:uex} \tag{30}$$$

which by insertion in the PDE leads to $$f(x,t)=2(1+t)c^2$$. This $$u_e$$ fulfills the boundary conditions $$u=0$$ and demands $$I(x)=x(L-x)$$ and $$V(x)={\frac{1}{2}}x(L-x)$$.

To realize that the chosen $$u_e$$ is also an exact solution of the discrete equations, we first remind ourselves that $$t_n=n\Delta t$$ so that

$$$\lbrack D_tD_t t^2\rbrack^n = \frac{t_{n+1}^2 - 2t_n^2 + t_{n-1}^2}{\Delta t^2} = (n+1)^2 -2n^2 + (n-1)^2 = 2, \label{_auto8} \tag{31}$$$
$$$\lbrack D_tD_t t\rbrack^n = \frac{t_{n+1} - 2t_n + t_{n-1}}{\Delta t^2} = \frac{((n+1) -2n + (n-1))\Delta t}{\Delta t^2} = 0 \label{_auto9} \tag{32}$$$

Hence,

$[D_tD_t u_e]^n_i = x_i(L-x_i)[D_tD_t (1+{\frac{1}{2}}t)]^n = x_i(L-x_i){\frac{1}{2}}[D_tD_t t]^n = 0$

Similarly, we get that

\begin{split} \begin{align*} \lbrack D_xD_x u_e\rbrack^n_i &= (1+{\frac{1}{2}}t_n)\lbrack D_xD_x (xL-x^2)\rbrack_i\\ & = (1+{\frac{1}{2}}t_n)\lbrack LD_xD_x x - D_xD_x x^2\rbrack_i \\ &= -2(1+{\frac{1}{2}}t_n) \end{align*} \end{split}

Now, $$f^n_i = 2(1+{\frac{1}{2}}t_n)c^2$$, which results in

$[D_tD_t u_e - c^2D_xD_xu_e - f]^n_i = 0 + c^2 2(1 + {\frac{1}{2}}t_{n}) + 2(1+{\frac{1}{2}}t_n)c^2 = 0$

Moreover, $$u_e(x_i,0)=I(x_i)$$, $$\partial u_e/\partial t = V(x_i)$$ at $$t=0$$, and $$u_e(x_0,t)=u_e(x_{N_x},0)=0$$. Also the modified scheme for the first time step is fulfilled by $$u_e(x_i,t_n)$$.

Therefore, the exact solution $$u_e(x,t)=x(L-x)(1+t/2)$$ of the PDE problem is also an exact solution of the discrete problem. This means that we know beforehand what numbers the numerical algorithm should produce. We can use this fact to check that the computed $$u^n_i$$ values from an implementation equals $$u_e(x_i,t_n)$$, within machine precision. This result is valid regardless of the mesh spacings $$\Delta x$$ and $$\Delta t$$! Nevertheless, there might be stability restrictions on $$\Delta x$$ and $$\Delta t$$, so the test can only be run for a mesh that is compatible with the stability criterion (which in the present case is $$C\leq 1$$, to be derived later).

Notice.

A product of quadratic or linear expressions in the various independent variables, as shown above, will often fulfill both the PDE problem and the discrete equations, and can therefore be very useful solutions for verifying implementations.

However, for 1D wave equations of the type $$u_{tt}=c^2u_{xx}$$ we shall see that there is always another much more powerful way of generating exact solutions (which consists in just setting $$C=1$$ (!), as shown in the section Analysis of the difference equations).


text:
# All Questions

3 views

### Multiple dynamic equations for plotting?

Say I have a bunch of equations between variables e.g. a=bcd, d=bc, e=bda and I want to do a plot between any two of them, is there an easier way than just rearranging by hand and then plugging into ...
24 views

### Use the trigonometric identities to expand and simplify if possible

Use the trigonometric identities to expand and simplify if possible. A, cos(270∘−θ) B, sin(270∘−θ) C, cos(270∘+θ) Here is my answers, it needs checking. 1. cos(270 - theta) = cos( 180 +90 - ...
58k views

### Where can I find examples of good Mathematica programming practice?

I consider myself a pretty good Mathematica programmer, but I'm always looking out for ways to either improve my way of doing things in Mathematica, or to see if there's something nifty that I haven't ...
81 views

### How to make hollow octahedron 3-compound?

I would like to remove facets from a graphic - like on the picture below. Link I tried to combine these two graphics: ...
24 views

### Interpolating a polynomial for with first and second derivative of points

I am trying to fit a polynomial for switch cutoff which is used in Molecular dynamics. I want a fifth order polynomial with six coefficient $c_0, ...,c_6$ between $r_1<r<r_2$ with these ...
3k views

### How to handle NDSolve::ndsz problem (singularity problem)

I have 2 second order differential equations (non-linear). The physics behind them is correct. I verified the equations many times. It is a solid pendulum with a mass-spring at the end of it. Now, ...
42 views

### Help with the plot f(x)=2ln(-0.5x+2)-1

I don't know how to write equation in the notebook. The graph doesn't appear when I tried.
35 views

### How can I graph the polynomial function? [on hold]

There is no explanation for the polynomial function. The equation is f(x)=-2x^2-4x-3.
125 views

### Coordinate free differential forms package

I am looking for a package in Mathematica which can handle differential forms in a coordinate free manner. I am aware of several packages which do differential forms, but it seems that for all of them ...
104 views

### Obtain polygons describing all intersections of many polygons

I have 6894 polygons describing zones in a state plane (i.e., cartesian) coordinate system. Most are not very complex, and I believe none have holes. A random example is Polygon[{{351633., ...
25 views

### How to keep first observation and delete other observations of each group in unbalanced Panel data

I have a huge panel data. I want to create a new data keeping only first observation of each group. How can I achieve it? For illustrative purpose; ...
35 views

### Copy files existing in different directories under the same directory

I have many files existing in different directories, I need to copy it under the same directory. Any suggestions?
20 views

### Mixture coefficients and Parameter estimation for moving Normal distributions

I am searching for a way to do the following steps below. I am stumped and I have stumped a few profs too. Any help would be appreciated. I would like to import data from, say, excel or a .txt file ...
96 views
+100

### Developing a function of two variables from given data

Cross listed from MSE: http://math.stackexchange.com/questions/1029399/developing-a-function-of-two-variables-from-given-data (I believe Mathematica SE is an appropriate place to ask this, as well.) ...
55 views

### How to tell Mathematica that certain variables are real/imaginary, integer-valued, etc

I'm trying to expedite some quantum mechanical calculations by running them through Mathematica. Expectation values and stuff like that. When I say, for example, ...
38 views

### The use of Evaluate in NDSolve and Plotting the answer [duplicate]

I am having a bit of hard time understanding the use of Evaluate when plotting the output of NDSolve. From the documentation, I enter: ...
27 views

### How do I include assignments inside this?

Noobie question, but if I have this: ...
39 views

### Default color for Filling in Mathematica 9

The default color scheme for plots in my Mathematica version is ColorData[1]. What about for the Filling option?
37 views

### DateListLogLogPlot

The command completion feature in Version 10.0.1.0 shows a command DateListLogLog but there is no documentation and doesn't work. Is this something that is going to be released in the near future and ...
56 views

### Coding While loops

I am interested in replicating the code in this question on MSE, whose contents I have copied for ease of reading: Pseudocode: Consider that 1 is the starting index of a list 1.  input natural ...
52 views

### Make plots of solutions persist

How can the code below be modified so that at allow plots of solutions to persist upon clicking in the slope field. For example, after three clicks, the plots of all three solutions should remain. In ...
117 views
+50

### How to use C language functions in Mathematica?

I'm trying to find out the way to use c-functions in Mathematica (I'm using 14 version) on OS X 10.10. For example I've got this lines of code: ...
170 views

### How to write Morse code into Plot Dashing?

How can we define a dashing pattern in a plot such that it reads arbitrary text in Morse code? For instance Graphics[{DotDashed, Line[{{0, 0}, {5, 0}}]}] ...
64 views

### PlotStyle and Legends

I know only three plot style ...
343 views

### Camel hump string split

Given System names of 2 or more characters: systemNames = Names["System*"] // Select[StringLength[#] > 1 &] ; Is ...
38 views

### What extrapolation method does Interpolation use by default?

As far as I can tell, Interpolation uses cubic spline interpolation by default - please correct me if I'm wrong. What does it use for extrapolation outside of the ...
52 views

### How to embed a filled Region in a Graphics?

Let for example: regn = RegionDifference[Rectangle[{-2, -2}, {2, 2}], Disk[{0, 0}, 1]]; The only documented way I found to draw a gereric region ...
21 views

92 views

### AstronomicalData not including units?

I like using mathematica to do my physics homework, and AstronomicalData is useful for that, but it doesn't seem to work properly with units. [the question is, at ...
138 views

### Combining solar system planets orbits and positions in Graphics3D and Show

I am building an animated model of the solar system, where I would like to display the OrbitPaths of the planets along with their ...
312 views

### How to get the distance of Jupiter from the sun?

I am aware of the command AstronomicalData["Jupiter", "Distance"] from which I can get the distance between the earth and Jupiter. Is there any way to use ...
103 views

This code: ...
34 views

### UnitConvert to Electronvolts converts to Electronvolts/Moles

The question probably is clear from the title. When I want to convert a kinetic energy, that contains the Molar Gas-Constant into Electronvolts I get eV/mol. ...
95 views

### Why does LyapunovSolve solve non-standard form?

LyapunovSolve[A,Q] solves the equation $A P + P A^T = Q$ for $P$ whereas the standard form (wikipedia,lecture notes, p.25, linear control systems) of the Lyapunov ...
### How to obtain Padé approximant of $\log(x)$?
Having read this answer on Math.SE, I wanted to try seeing how Padé approximations converge to $\log(x)$. But my first attempt after reading the documentation failed: ...

format: plain


text:
Biostar Beta. Not for public use.
common sequence between multi fasta file
0
Entering edit mode
17 months ago
Sam • 110

Dear All

how I can retrieve common seqs between multiple fasta (12 ) files?

Thanks

bash • 732 views
1
Entering edit mode

I already obtain common seqs with this command , is it a correct way ?

comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(comm -12 <(sort A1_novel_mirnas.fa) <(sort A2_novel_mirnas.fa)) <(sort B1_novel_mirnas.fa)) <(sort B2_novel_mirnas.fa)) <( sort C1_novel_mirnas.fa)) <(sort C2_novel_mirnas.fa)) <(sort D1_novel_mirnas.fa)) <(sort D2_novel_mirnas.fa)) <(sort E1_novel_mirnas.fa)) <(sort E2_novel_mirnas.fa)) <(sort F1_novel_mirnas.fa)) <(sort F2_novel_mirnas.fa)

0
Entering edit mode

So you have RNA sequence. Can you try replacing U with T and running @5heikki's script again?

0
Entering edit mode

Hello Sam,

what do you mean by "common seqs". Please describe more clearly what you have and what your goal is.

fin swimmer

0
Entering edit mode

Hi

Common mean sequence which is available in all fasta file. I want to find common sequences between all 12 fasta files.

0
Entering edit mode

Are there many sequences in each file or just one per file?

0
Entering edit mode

1
Entering edit mode

Doing quick blat searches is one option that should be pretty fast. Other option may be CD-HIT to do pair-wise comparisons.

0
Entering edit mode

The unix commandline utility comm will be useful here, so long as you can sort your files first. (Assuming ‘common’ in this question means the sequences or headers are identical)

0
Entering edit mode
1. make a list of the unique fasta sequences

cat *.mfa | grep -v ">" | sort -u > unique_sequences

2. loop each line in all .mfa file

cat unique_sequences | while read line; do grep -w $line *.mfa; done 3. you could easily go further and check which sequence is present 12 times ADD REPLYlink 1 Entering edit mode 18 months ago 5heikki 8.4k Finland This assumes that 1) blast is installed and that 2) the script is run in the dir where .fasta named files reside #!/bin/bash FILES=($(find . -maxdepth 1 -type f -name "*.fasta"))
for((i=0;i<${#FILES[@]}-1;i++)); do blastn -query${FILES[$i]} -subject${FILES[$i+1]} -outfmt '6 qseqid qlen length slen qseq' \ | awk 'BEGIN{OFS=FS="\t"}{if($2==$3 &&$2==$4){print ">"$1"\n"$5}}' > tmp.fasta cp tmp.fasta${FILES[\$i+1]}
done

So we blast first file against second file and put into tmp.fasta only the identical sequences between the files. Then we replace the second file with tmp.fasta. Then we blast that against the third file. Etc. Files will be deleted. If there are common sequences, in the end they will be in the tmp.fasta file. Here we assume that if query length = alignment length = subject length the sequences are the same. Of course, this might not be the case. I leave it for you how to make it more precise blastn -help

0
Entering edit mode

Did you mean to use blast2seq? I don't see blast indexes being made anywhere. Edit: That is the correct syntax for blast2seq (not having used it on the command line made me check again). @sam: If you are getting an error then post the exact error or explain what "is not working".

0
Entering edit mode

unfortunately it's not working and as Genomax mentioned for blastn we need index , which is missing

0
Entering edit mode

If your file names don't end in .fasta then you need to change the code above to account for that.

0
Entering edit mode

Warning: [blastn] Query is Empty! and it convert all fasta file to empty file!

0
Entering edit mode

I just tested @5heikki's script and it worked flawlessly.

Do all your files have consistent names/format?

0
Entering edit mode

If the sequences are really short add -task blastn-short to the command. Also, is it possible that you don't have any "universal" seqs?

0
Entering edit mode

the length of seqs are about 20 bp , how ever I tested with blastn-short and it produced a huge temp file over than 600 Mb from just 600 seqs (total amount of seqs from all files is about 600 ) . when I lunched it with just 4 lib it was fine but for all 12 libs it's not work correctly

Check from blastn -help output how to add a field to outfmt for sequence similarity and filter also based on that in the awk command. In case the "temp" file is still huge, you must have a lot of redundancy in your input.


text:
5

# NOTE: This IS 3 mutti-part question Once an- Block Ahas mass of 65 kg; and answer IS submitted You will be unable to return to this part block Bhas are Us mass 0f 1...

## Question

###### NOTE: This IS 3 mutti-part question Once an- Block Ahas mass of 65 kg; and answer IS submitted You will be unable to return to this part block Bhas are Us mass 0f 13kg 0.20 and Uk= 015. The coefficients of friction between all sunaces of contactdetermlne the acceleration ol block(You must provlde an answer before moving on to the next part )The acceleratlon cf block B Isms2 + 25

NOTE: This IS 3 mutti-part question Once an- Block Ahas mass of 65 kg; and answer IS submitted You will be unable to return to this part block Bhas are Us mass 0f 13kg 0.20 and Uk= 015. The coefficients of friction between all sunaces of contact determlne the acceleration ol block (You must provlde an answer before moving on to the next part ) The acceleratlon cf block B Is ms2 + 25

#### Similar Solved Questions

##### Suppose that the mechanism of hydrolysi8 of tert-butyl chloride proceeded by & different mechaniam 8 that the predicted rate law was: JRCH KRRCIJHHzO]Fron this experiment, would you be able to distinguish between this mechaniam (which gives an order reaction of reaction 2) and the mechanism proposed above (which gives firat order rate law)? Explain your angwer:
Suppose that the mechanism of hydrolysi8 of tert-butyl chloride proceeded by & different mechaniam 8 that the predicted rate law was: JRCH KRRCIJHHzO] Fron this experiment, would you be able to distinguish between this mechaniam (which gives an order reaction of reaction 2) and the mechanism pro...
##### At least one of the answers above is NOT correctpoint) Find the average rate of change of f(z) 312 + between each of the pairs of points below:(a) Between (2,13) and (5,76)answer5607(b) Between (b, 9) and (9, w)answer3(g+w)(g-w)(b-q)(c) Between (â‚¬, f(z)) ad (2 + h, f(z + h))answer6x+3h
At least one of the answers above is NOT correct point) Find the average rate of change of f(z) 312 + between each of the pairs of points below: (a) Between (2,13) and (5,76) answer 5607 (b) Between (b, 9) and (9, w) answer 3(g+w)(g-w)(b-q) (c) Between (â‚¬, f(z)) ad (2 + h, f(z + h)) answer 6x+...
##### Rittu Kalia Assignment Math[HWeh Work Task5 due (L/22/2uzu at [:OOpm PSTufv-mathllI-Schlitt1. (1 pint) Evaluate the limit7. (1 point) Evaluate the limitr_81 lim "950- |38(lf the limit does not exist enter DNE" . ! Take un algebraie approch" ( this (so #nol* using YuM calculator) USIDg Elementary Function Theorem in elass Practice justilying evet; step, DeM though the problem just asks for numerical answer: Limil 8. ( [ point) Reduce the expression and then evaluate the limitTake
Rittu Kalia Assignment Math[HWeh Work Task5 due (L/22/2uzu at [:OOpm PST ufv-mathllI-Schlitt 1. (1 pint) Evaluate the limit 7. (1 point) Evaluate the limit r_81 lim "950- |38 (lf the limit does not exist enter DNE" . ! Take un algebraie approch" ( this (so #nol* using YuM calculator...
##### Consider the AR(I) process given by0.6yt-1where et is a white noise process with Elet) 0 and Vlet) Ele 2) = 1.5 for all t Substitute Yt-1 with (0.6yt--2 et-1), then Yt-2 with (0.6yt-3 et-2), and so forth: Assume that the starting value of {Yt}, which is denoted Yo, is zero(1) What is the effect of a random shock (of an increase by one unit) that occurred three periods ago (et--3) on Yt?(2) What is Elyt)?(3) If Cov(yt-1, er) = 0and V(yt-1) = Vly ) = 02 for all t, 02 is(4) What is Cov(yt et)?(5) W
Consider the AR(I) process given by 0.6yt-1 where et is a white noise process with Elet) 0 and Vlet) Ele 2) = 1.5 for all t Substitute Yt-1 with (0.6yt--2 et-1), then Yt-2 with (0.6yt-3 et-2), and so forth: Assume that the starting value of {Yt}, which is denoted Yo, is zero (1) What is the effect o...
##### Calculate the total error by finding the model equation in the form y = 1 Ae" that fits the points y) with the least error in the table below:1.52.55.50.380.430.680.790.900.940.960.970.980.99
Calculate the total error by finding the model equation in the form y = 1 Ae" that fits the points y) with the least error in the table below: 1.5 2.5 5.5 0.38 0.43 0.68 0.79 0.90 0.94 0.96 0.97 0.98 0.99...
##### CompleteReviewPart ADetermine the position (or positions) (Figure 1) where an electron would experience no net force. there is more than one answer; enter them in ascending order separated by commals): If there is no answer; enterFigureof 1eubrPrevious AnswersCorrectProvide FeedbackNext
Complete Review Part A Determine the position (or positions) (Figure 1) where an electron would experience no net force. there is more than one answer; enter them in ascending order separated by commals): If there is no answer; enter Figure of 1 eubr Previous Answers Correct Provide Feedback Next...
##### Which carboxylic acid contains three hydroxyl groups: salicylic acid Ballic acid c. citric acid pyruvic acid
Which carboxylic acid contains three hydroxyl groups: salicylic acid Ballic acid c. citric acid pyruvic acid...
##### QUESTION 6 Fill in the blanks for the statoments below s0 that they may read as truestoreoisomer Ihat has mirror image whichnon-superimposablemixture af 2 onantiomers in a 1 1 ratio_ The conformation tho most stable structuro of cyclohexane Isorers are described Isomers that have different connactivity but the same molecular formula_ chiral carbon also known as Muy be used measure (he rotation plane-palarized Iight (PPL) substance Ihat rotatos PPL clockwise direction said[7x2-14]ORGANIC CH
QUESTION 6 Fill in the blanks for the statoments below s0 that they may read as true storeoisomer Ihat has mirror image which non- superimposable mixture af 2 onantiomers in a 1 1 ratio_ The conformation tho most stable structuro of cyclohexane Isorers are described Isomers that have different c...
##### Inaider the icllowing s0l cf dependunt and indepondenl vanables Complato parts . Ihrough belJolormina nppraprale crilical value(a) Ix u-0.10 ,onvca valuo(a) ialato) Round Iwn docimal placcunoodod coiner wpnrule uneters 09 F neded ) Srlo ira conclubion lot to toel [ "igndlicancn Iur *Io nu] hypothonin ThoroWdeucconclde Ihat 4ntorolatonahip DottconUayInkorprot Ita p-waluo Ior ench Indopor dont vnrablo Intexptat Ina Henkn hypobiue? "at cort paponceng Io * wnd 4-0,10Sinco Uhut Ihu nud ny
Inaider the icllowing s0l cf dependunt and indepondenl vanables Complato parts . Ihrough bel Jolormina nppraprale crilical value(a) Ix u-0.10 , onvca valuo(a) ialato) Round Iwn docimal placcu noodod coiner wpnrule uneters 09 F neded ) Srlo ira conclubion lot to toel [ "igndlicancn Iur * Io nu]...
##### Face of the Tinees the front face house height ol the truss. 1 asnou W The M W 5 Ju post rectangular 2 region. than the height of the truss above; Round - The eruesgth nearest If the ehetze area Teclaucula( of the eegiont !The 1 length Snougpai degion region is
face of the Tinees the front face house height ol the truss. 1 asnou W The M W 5 Ju post rectangular 2 region. than the height of the truss above; Round - The eruesgth nearest If the ehetze area Teclaucula( of the eegiont ! The 1 length Snougpai degion region is...
##### Write each expression as the sine, cosine, or tangent of a double angle. Then find the exact value of the expression. $$2 \cos ^{2} \frac{\pi}{8}-1$$
Write each expression as the sine, cosine, or tangent of a double angle. Then find the exact value of the expression. $$2 \cos ^{2} \frac{\pi}{8}-1$$...
##### Waves on a Stick. A flexible stick 2.0 $\mathrm{m}$ long is not fixed in any way and is free to vibrate. Make clear drawings of this stick vibrating in its first three harmonics, and then use your drawings to find the wavelengths of each of these harmonics. (Hint: Should the ends be nodes or antinodes?)
Waves on a Stick. A flexible stick 2.0 $\mathrm{m}$ long is not fixed in any way and is free to vibrate. Make clear drawings of this stick vibrating in its first three harmonics, and then use your drawings to find the wavelengths of each of these harmonics. (Hint: Should the ends be nodes or antinod...
##### Suppose a mutation occurred that inactivated Rab5-GAP. As aresult of this mutation, which of the following would you expect tosee an increase in the concentration of:Rab5-GTPKDELmitochondrial DNAmannose-6-phosphate (M6P)Rab5-GDP
Suppose a mutation occurred that inactivated Rab5-GAP. As a result of this mutation, which of the following would you expect to see an increase in the concentration of: Rab5-GTP KDEL mitochondrial DNA mannose-6-phosphate (M6P) Rab5-GDP...
##### Describe what is meant by Big Crunch, Big Chill,and Big Rip.
Describe what is meant by Big Crunch, Big Chill, and Big Rip....


text:
Volume 282 - 38th International Conference on High Energy Physics (ICHEP2016) - Higgs Physics
Measurement of cross sections and couplings of the SM Higgs boson in the ZZ decay channel using the ATLAS detector
L. Aperio Bella* on behalf of the ATLAS collaboration
*corresponding author
Full text: pdf
Pre-published on: February 06, 2017
Published on: April 19, 2017
Abstract
The Higgs boson decay to four leptons, provides good
sensitivity for the measurement of its properties due to its high signal-to-background
ratio (S/B), which is better than two for each of the four final states: $4\mu$,
$2e2\mu$, $2\mu2e$ and $4e$, where the
first lepton pair is defined to be the one with the dilepton invariant mass closest to the $Z$ boson mass.
These proceedings summarise several preliminary results on the properties of the Higgs boson decaying to four leptons final state
using $14.8$~fb$^{-1}$ of integrated luminosity collected in $pp$ collisions by the ATLAS detector at the LHC
at a centre-of-mass energy of 13 TeV. The results include measurements of the fiducial cross sections, the different
production mode cross sections and a test of beyond the Standard Model interactions using
an effective Lagrangian approach.
DOI: https://doi.org/10.22323/1.282.0386
How to cite

Metadata are provided both in "article" format (very similar to INSPIRE) as this helps creating very compact bibliographies which can be beneficial to authors and readers, and in "proceeding" format which is more detailed and complete.

Open Access


text:
# Difference between revisions of "Math 120"

Hi Matthew, let's see if this test works: Here's some random latex code: $\displaystyle f(x)=\left(x^{3}-2 x\right)^{3}, \hspace{2mm} f^{\prime \prime}(1)$


text:
# ISO: cross product

ISO norms specifies a standard symbol for cross product? I found in many books of mathematics the

\times

operator, while in physics it seems widely used the

 \wedge

(or \land or another similar symbol?) operator.

Which is the correct one to ISO norms?

-
As far as I know, both conventions are used in different communities. –  Corentin Oct 15 '12 at 19:59
@Corentin Of course. This was a "personal" statistics and it is meaningless for the question. ; ) (It was just not to put the two commands only.) –  R. M. Oct 15 '12 at 20:03

The standard ISO 80000-2, “Mathematical signs and symbols to be used in the natural sciences and technology”, specifies that the symbol for vector product (commonly known as cross product, too) is “×” MULTIPLICATION SIGN U+00D7. This corresponds to \times in LaTeX.

The “wedge” symbol “∧” denotes “logical and” by that standard and in common practice in formal logic. Using it for other purposes is thus nonstandard and potentially misleading.

-
Well, it's quite standard to use \wedge when talking about the wedge (exterior) product. –  Ipsen Oct 15 '12 at 20:58
The “wedge” symbol “∧” is used as the “logical and”, but it is also used by many authors as the "wedge or exterior product" of two vectors. –  Américo Tavares Oct 15 '12 at 21:00
\wedge is also used in Topology, to denote the quotient space $X\wedge Y=X\times Y / X\vee Y$, called the smash product of $X$ and $Y$. –  Sigur Oct 16 '12 at 0:41
@Sigur Mathematicians are not subject to that ISO regulation, and no regulation would be able to tame them. The \wedge symbol is used for very different operations, for example it's widely used as a connective in formal logic ("and") or to denote the infimum in lattices (e.g., Boolean algebras). The ISO standard is for what its title specifies: math symbols to be used in natural sciences and technology. –  egreg Oct 16 '12 at 6:52
You can use any symbol you want for any mathematical action you want, as long as you define it properly. In "tropical algebras" you say "multiplication is addition and addition is minimization". Btw, \wedge and \vee are used not only in logic, but in a whole theory of lattices, as sort-of generalization of and and or. And ISO standard cannot stop physicists neither, because they need new and new things that no standard can cover. –  tohecz Oct 16 '12 at 7:29
show 1 more comment

u × v means the cross product of vectors u and v.

u ∧ v means the wedge product of any multivectors u and v. In three dimensional Euclidean space the wedge product and the cross product of two vectors are each other's Hodge dual.

So, × (\times) is the cross product and ∧ (\wedge) is the wedge product.

-
I didn't know there were any difference... Why they are used as the same thing in R^3? Which is the correct one to use in physics? –  R. M. Oct 15 '12 at 20:38
@R.M. Well, I would ask this question in physics.SE. In the meantime I saw that in the German version of the List of mathematical symbols the two symbols are called cross product (Kreuzprodukt) –  Américo Tavares Oct 15 '12 at 20:48
Because they basically mean the same thing on tensors of dimension one (=vectors). Their difference is in higher dimension tensors. –  tohecz Oct 16 '12 at 7:23


text:
A library that integrates huggingface transformers with version 2 of the fastai framework

Install

You can now pip install blurr via pip install ohmeow-blurr

Or, even better as this library is under very active development, create an editable install like this:

git clone https://github.com/ohmeow/blurr.git
cd blurr
pip install -e ".[dev]"

How to use

The initial release includes everything you need for sequence classification and question answering tasks. Support for token classification and summarization are incoming. Please check the documentation for more thorough examples of how to use this package.

The following two packages need to be installed for blurr to work:

1. fastai2 (see http://docs.fast.ai/ for installation instructions)
2. huggingface transformers (see https://huggingface.co/transformers/installation.html for details)

Imports

import torch
from transformers import *
from fastai.text.all import *

from blurr.data.all import *
from blurr.modeling.all import *

[nltk_data] Downloading package wordnet to /home/wgilliam/nltk_data...
[nltk_data]   Package wordnet is already up-to-date!

path = untar_data(URLs.IMDB_SAMPLE)

model_path = Path('models')

task = HF_TASKS_AUTO.SequenceClassification

pretrained_model_name = "bert-base-uncased"

blocks = (HF_TextBlock(hf_arch, hf_config, hf_tokenizer, hf_model), CategoryBlock)

dls.show_batch(dataloaders=dls, max_n=2)

text category
0 raising victor vargas : a review < br / > < br / > you know, raising victor vargas is like sticking your hands into a big, steaming bowl of oatmeal. it's warm and gooey, but you're not sure if it feels right. try as i might, no matter how warm and gooey raising victor vargas became i was always aware that something didn't quite feel right. victor vargas suffers from a certain overconfidence on the director's part. apparently, the director thought that the ethnic backdrop of a latino family on the lower east side, and an idyllic storyline would make the film critic proof. he was right, but it didn't fool me. raising victor vargas is the story about a seventeen - year old boy called, you guessed it, victor vargas ( victor rasuk ) who lives his teenage years chasing more skirt than the rolling stones could do in all the years they've toured. the movie starts off in  ugly fat'donna's bedroom where victor is sure to seduce her, but a cry from outside disrupts his plans when his best - friend harold ( kevin rivera ) comes - a - looking for him. caught in the attempt by harold and his sister, victor vargas runs off for damage control. yet even with the embarrassing implication that he's been boffing the homeliest girl in the neighborhood, nothing dissuades young victor from going off on the hunt for more fresh meat. on a hot, new york city day they make way to the local public swimming pool where victor's eyes catch a glimpse of the lovely young nymph judy ( judy marte ), who's not just pretty, but a strong and independent too. the relationship that develops between victor and judy becomes the focus of the film. the story also focuses on victor's family that is comprised of his grandmother or abuelita ( altagracia guzman ), his brother nino ( also played by real life brother to victor, silvestre rasuk ) and his sister vicky ( krystal rodriguez ). the action follows victor between scenes with judy and scenes with his family. victor tries to cope with being an oversexed pimp - daddy, his feelings for judy and his grandmother's conservative catholic upbringing. < br / > < br / > the problems that arise from raising victor vargas are a few, but glaring errors. throughout the film you get to know certain characters like vicky, nino, grandma, judy and even negative
1 i really wanted to love this show. i truly, honestly did. < br / > < br / > for the first time, gay viewers get their own version of the " the bachelor ". with the help of his obligatory " hag " andra, james, a good looking, well - to - do thirty - something has the chance of love with 15 suitors ( or " mates " as they are referred to in the show ). the only problem is half of them are straight and james doesn't know this. if james picks a gay one, they get a trip to new zealand, and if he picks a straight one, straight guy gets \$ 25, 000. how can this not be fun?! take my hand, lets stroll : < br / > < br / > the most glaring problem with this show is the bachelor himself. james is your typical young and successful gay guy with a nice smile and body, the one you'd probably give two glances towards at your local bar before grazing for greener pastures. why they chose to cast james as the leading man is beyond me. god knows there's so many other hotter and vivacious homosexual men out there dying to be on tv. < br / > < br / > aside from his rather average physical appearance, james is about as interesting and exciting as a piece of chalk. even as such, he has this arrogant, smugly condescending aura about him. however, if james were standing up against a blank, white wall he'd meld right into in it. i honestly can't recall a single interesting or noteworthy thing james said during the course of the show. he is that boring and forgettable. in fact, one of the mates flat out advised him he wasn't feeling a connection. i thought that was the best part of the show. also, james speaks with an excruciatingly annoying lilt. sound feminine or sound masculine, but don't * * * * ing segue tones in the middle of sentences... so painful to sit through. i hated him so much all throughout the show i kept thinking, " please choose a straight guy and humiliate yourself and your unfortunate looking hag " < br / > < br / > then we have the suitors. a remarkably bland bunch of men who don't seem to care either way what is happening. equally vapid, they seem to be indistinguishable from one guy to the next except, negative

... and 🚂

model = HF_BaseModelWrapper(hf_model)

learn = Learner(dls,
model,
loss_func=CrossEntropyLossFlat(),
metrics=[accuracy],
cbs=[HF_BaseModelCallback],
splitter=hf_splitter)

learn.create_opt()
learn.freeze()

learn.fit_one_cycle(3, lr_max=1e-3)

epoch train_loss valid_loss accuracy time
0 0.649641 0.768956 0.585000 00:22
1 0.385529 0.295318 0.895000 00:22
2 0.318808 0.284280 0.885000 00:22
learn.show_results(learner=learn, max_n=2)

text category target

12/31/2020

The "Goodbye 2020" release with lots of goodies for blurr users:

• Updated the Seq2Seq models to use some of the latest huggingface bits like tokenizer.prepare_seq2seq_batch.
• Separated out the Seq2Seq and Token Classification metrics into metrics-specific callbacks for a better separation of concerns. As a best practice, you should now only use them as fit_one_cycle, etc.. callbacks rather than attach them to your Learner.
• NEW: Translation are now available in blurr, joining causal language modeling and summarization in our core Seq2Seq stack
• NEW: Integration of huggingface's Seq2Seq metrics (rouge, bertscore, meteor, bleu, and sacrebleu). Plenty of info on how to set this up in the docs.
• NEW: Added default_text_gen_kwargs, a method that given a huggingface config, model, and task (optional), will return the default/recommended kwargs for any text generation models.
• A lot of code cleanup (e.g., refactored naming and removal of redundant code into classes/methods)
• More model support and more tests across the board! Check out the docs for more info
• Misc. validation improvements and bug fixes.

As I'm sure there is plenty I can do to make this library better, please don't hesitate to join in and help the effort by submitting PRs, pointing out problems with my code, or letting me know what and how I can improve things generally. Some models, like mbart and mt5 for example, aren't giving good results and I'd love to get any and all feedback from the community on how to resolve such issues ... so hit me up, I promise I won't bit :)

12/20/2020

12/12/2020

• Updated to work with the latest version of fast.ai (2.1.8) and huggingface transformers >= 4.x
• Fixed Learner.blurr_summary to work with fast.ai >= 2.1.8
• Fixed inclusion of add_prefix_space in tokenizer BLURR_MODEL_HELPER
• Fixed token classification show_results for tokenizers that add a prefix space
• Notebooks run with environment variable "TOKENIZERS_PARALLELISM=false" to avoid fast tokenizer warnings
• Updated docs

11/12/2020

• Updated documentation
• Updated model callbacks to support mixed precision training regardless of whether you are calculating the loss yourself or letting huggingface do it for you.

11/10/2020

• Major update just about everywhere to facilitate a breaking change in fastai's treatment of before_batch transforms.
• Reorganized code as I being to work on LM and other text2text tasks
• Misc. fixes

10/08/2020

• Updated all models to use ModelOutput classes instead of traditional tuples. ModelOutput attributes are assigned to the appropriate fastai bits like Learner.pred and Learner.loss and anything else you've requested the huggingface model to return is available via the Learner.blurr_model_outputs dictionary (see next two bullet items)
• Added ability to grab attentions and hidden state from Learner. You can get at them via Learner.blurr_model_outputs dictionary if you tell HF_BaseModelWrapper to provide them.
• Added model_kwargs to HF_BaseModelWrapper should you need to request a huggingface model to return something specific to it's type. These outputs will be available via the Learner.blurr_model_outputs dictionary as well.

09/16/2020

• Major overhaul to do everything at batch time (including tokenization/numericalization). If this backfires, I'll roll everything back but as of now, I think this approach not only meshes better with how huggingface tokenization works and reduce RAM utilization for big datasets, but also opens up opportunities for incorporating augmentation, building adversarial models, etc.... Thoughts?
• Added tests for summarization bits
• New change may require some small modifications (see docs or ask on issues thread if you have problems you can't fiture out). I'm NOT doing a release until pypi until folks have a chance to work with the latest.

09/07/2020

• Updated summarization to support BART, T5, and Pegasus

08/20/2020

• Updated everything to work latest version of fastai (tested against 2.0.0)
• Added batch-time padding, so that by default now, HF_TokenizerTransform doesn't add any padding tokens and all huggingface inputs are padded simply to the max sequence length in each batch rather than to the max length (passed in and/or acceptable to the model). This should create efficiencies across the board, from memory consumption to GPU utilization. The old tried and true method of padding during tokenization requires you to pass in padding='max_length to HF_TextBlock.
• Removed code to remove fastai2 @patched summary methods which had previously conflicted with a couple of the huggingface transformers

08/13/2020

• Updated everything to work latest transformers and fastai
• Reorganized code to bring it more inline with how huggingface separates out their "tasks".

07/06/2020

• Updated everything to work huggingface>=3.02
• Changed a lot of the internals to make everything more efficient and performant along with the latest version of huggingface ... meaning, I have broken things for folks using previous versions of blurr :).

06/27/2020

• Simplified the BLURR_MODEL_HELPER.get_hf_objects method to support a wide range of options in terms of building the necessary huggingface objects (architecture, config, tokenizer, and model). Also added cache_dir for saving pre-trained objects in a custom directory.
• Misc. renaming and cleanup that may break existing code (please see the docs/source if things blow up)
• Added missing required libraries to requirements.txt (e.g., nlp)

05/23/2020

• Initial support for text generation (e.g., summarization, conversational agents) models now included. Only tested with BART so if you try it with other models before I do, lmk what works ... and what doesn't

05/17/2020

• Major code restructuring to make it easier to build out the library.
• HF_TokenizerTransform replaces HF_Tokenizer, handling the tokenization and numericalization in one place. DataBlock code has been dramatically simplified.
• Tokenization correctly handles huggingface tokenizers that require add_prefix_space=True.
• HF_BaseModelCallback and HF_BaseModelCallback are required and work together in order to allow developers to tie into any callback friendly event exposed by fastai2 and also pass in named arguments to the huggingface models.
• show_batch and show_results have been updated for Question/Answer and Token Classification models to represent the data and results in a more easily intepretable manner than the defaults.

05/06/2020

• Initial support for Token classification (e.g., NER) models now included
• Extended fastai's Learner object with a predict_tokens method used specifically in token classification
• HF_BaseModelCallback can be used (or extended) instead of the model wrapper to ensure your inputs into the huggingface model is correct (recommended). See docs for examples (and thanks to fastai's Sylvain for the suggestion!)
• HF_Tokenizer can work with strings or a string representation of a list (the later helpful for token classification tasks)
• show_batch and show_results` methods have been updated to allow better control on how huggingface tokenized data is represented in those methods

⭐ Props

A word of gratitude to the following individuals, repos, and articles upon which much of this work is inspired from:


text:
src/HOL/Real/HahnBanach/Linearform.thy
 author haftmann Fri Jun 17 16:12:49 2005 +0200 (2005-06-17) changeset 16417 9bc16273c2d4 parent 14254 342634f38451 child 23378 1d138d6bb461 permissions -rw-r--r--
migrated theory headers to new format
     1 (*  Title:      HOL/Real/HahnBanach/Linearform.thy

     2     ID:         $Id$

     3     Author:     Gertrud Bauer, TU Munich

     4 *)

     5

     6 header {* Linearforms *}

     7

     8 theory Linearform imports VectorSpace begin

     9

    10 text {*

    11   A \emph{linear form} is a function on a vector space into the reals

    12   that is additive and multiplicative.

    13 *}

    14

    15 locale linearform = var V + var f +

    16   assumes add [iff]: "x \<in> V \<Longrightarrow> y \<in> V \<Longrightarrow> f (x + y) = f x + f y"

    17     and mult [iff]: "x \<in> V \<Longrightarrow> f (a \<cdot> x) = a * f x"

    18

    19 declare linearform.intro [intro?]

    20

    21 lemma (in linearform) neg [iff]:

    22   includes vectorspace

    23   shows "x \<in> V \<Longrightarrow> f (- x) = - f x"

    24 proof -

    25   assume x: "x \<in> V"

    26   hence "f (- x) = f ((- 1) \<cdot> x)" by (simp add: negate_eq1)

    27   also from x have "... = (- 1) * (f x)" by (rule mult)

    28   also from x have "... = - (f x)" by simp

    29   finally show ?thesis .

    30 qed

    31

    32 lemma (in linearform) diff [iff]:

    33   includes vectorspace

    34   shows "x \<in> V \<Longrightarrow> y \<in> V \<Longrightarrow> f (x - y) = f x - f y"

    35 proof -

    36   assume x: "x \<in> V" and y: "y \<in> V"

    37   hence "x - y = x + - y" by (rule diff_eq1)

    38   also have "f ... = f x + f (- y)" by (rule add) (simp_all add: x y)

    39   also from _ y have "f (- y) = - f y" by (rule neg)

    40   finally show ?thesis by simp

    41 qed

    42

    43 text {* Every linear form yields @{text 0} for the @{text 0} vector. *}

    44

    45 lemma (in linearform) zero [iff]:

    46   includes vectorspace

    47   shows "f 0 = 0"

    48 proof -

    49   have "f 0 = f (0 - 0)" by simp

    50   also have "\<dots> = f 0 - f 0" by (rule diff) simp_all

    51   also have "\<dots> = 0" by simp

    52   finally show ?thesis .

    53 qed

    54

    55 end


text:
If you're seeing this message, it means we're having trouble loading external resources on our website.

If you're behind a web filter, please make sure that the domains *.kastatic.org and *.kasandbox.org are unblocked.

Main content

# More examples of special products

CCSS.Math:

## Video transcript

And now I want to do a bunch of examples dealing with probably the two most typical types of polynomial multiplication that you'll see, definitely, in algebra. And the first is just squaring a binomial. So if I have x plus 9 squared, I know that your temptation is going to say, oh, isn't that x squared plus 9 squared? And I'll say, no, it isn't. You have to resist every temptation on the planet to do this. It is not x squared plus 9 squared. Remember, x plus 9 squared, this is equal to x plus 9, times x plus 9. This is a multiplication of this binomial times itself. You always need to remember that. It's very tempting to think that it's just x squared plus 9 squared, but no, you have to expand it out. And now that we've expanded it out, we can use some of the skills we learned in the last video to actually multiply it. And just to show you that we can do it in the way that we multiplied the trinomial last time, let's multiply x plus 9, times x plus a magenta 9. And I'm doing it this way just to show you when I'm multiplying by this 9 versus this x. But let's just do it. So we go 9 times 9 is 81. Put it in the constants' place. 9 times x is 9x. Then we have-- go switch to this x term-- we have a yellow x. x times 9x is 9x. Put it in the first degree space. x times x is x squared. And then we add everything up. And we get x squared plus 18x plus 81. So this is equal to x squared plus 18x plus 81. Now you might see a little bit of a pattern here, and I'll actually make the pattern explicit in a second. But when you square a binomial, what happened? You have x squared. You have this x times this x, gives you x squared. You have the 9 times the 9, which is 81. And then you have this term here which is 18x. How did we get that 18x? Well, we multiplied this x times 9 to get 9x, and then we multiplied this 9 times x to get another 9x. And then we added the two right here to get 18x. So in general, whenever you have a squared binomial-- let me do it this way. I'll do it in very general terms. Let's say we have a plus b squared. Let me multiply it this way again, just to give you the hang of it. This is equal to a plus b, times a plus-- I'll do a green b right there. So we have to b times b is b squared. Let's just assume that this is a constant term. I'll put it in the b squared right there. I'm assuming this is constant. So this would be a constant, this would be analogous to our 81. a is a variable that we-- actually let me change that up even better. Let me make this into x plus b squared, and we're assuming b is a constant. So it would be x plus b, times x plus a green b, right there. So assuming b's a constant, b times b is b squared. b times x is bx. And then we'll do the magenta x. x times b is bx. And then x times x is x squared. So when you add everything, you're left with x squared plus 2bx, plus b squared. So what you see is, the end product, what you have when you have x plus b squared, is x squared, plus 2 times the product of x and b, plus b squared. So given that pattern, let's do a bunch more of these. And I'm going to do it the fast way. So 3x minus 7 squared. Let's just remember what I told you. Just don't remember it, in the back of your mind, you should know why it makes sense. If I were to multiply this out, do the distributive property twice, you know you'll get the same answer. So this is going to be equal to 3x squared, plus 2 times 3x, times negative 7. Right? We know that it's 2 times each the product of these terms, plus negative 7 squared. And if we use our product rules here, 3x squared is the same thing as 9x squared. This right here, you're going to have a 2 times a 3, which is 6, times a negative 7, which is negative 42x. And then a negative 7 squared is plus 49. That was the fast way. And just to make sure that I'm not doing something bizarre, let me do it the slow way for you. 3x minus 7, times 3x minus 7. Negative 7 times negative 7 is positive 49. Negative 7 times 3x is negative 21x. 3x times negative 7 is negative 21x. 3x times 3x is 9 x squared. Scroll to the left a little bit. Add everything. You're left with 9x squared, minus 42x, plus 49. So we did indeed get the same answer. Let's do one more, and we'll do it the fast way. So if we have 8x minus 3-- actually, let me do one which has more variables in it. Let's say we had 4x squared plus y squared, and we wanted to square that. Well, same idea. This is going to be equal to this term squared, 4x squared, squared, plus 2 times the product of both terms, 2 times 4x squared times y squared, plus y squared, this term, squared. And what's this going to be equal to? This is going to be equal to 16-- right, 4 squared is 16-- x squared, squared, that's 2 times 2, so it's x to the fourth power. And then plus, 2 times 4 times 1, that's 8x squared y squared. And then y squared, squared, is y to the fourth. Now, we've been dealing with squaring a binomial. The next example I want to show you is when I take the product of a sum and a difference. And this one actually comes out pretty neat. So I'm going to do a very general one for you. Let's just do a plus b, times a minus b. So what's this going to be equal to? This is going to be equal to a times a-- let me make these actually in different colors-- so a minus b, just like that. So it's going to be this green a times this magenta a, a times a, plus, or maybe I should say minus, the green a times this b. I got the minus from right there. And then we're going to have the green b, so plus the green b times the magenta a. I'm just multiplying every term by every term. And then finally minus the green b-- that's where the minus is coming from-- minus the green b times the magenta b. And what is this going to be equal to? This is going to be equal to a squared, and then this is minus ab. This could be rewritten as plus ab, and then we have minus b squared. These right here cancel out, minus ab plus ab, so you're just left with a squared minus b squared. Which is a really neat result because it really simplifies things. So let's use that notion to do some multiplication. So if we say 2x minus 1, times 2x plus 1. Well, these are the same thing. The 2x plus 1, you could view this as, if you like, a plus b, and the 2x minus 1, you can view it as a minus b, where this is a, and that b is 1. This is b. That is a. Just using this pattern that we figured out just now. So what is this going to be equal to? It's going to be a squared, it's going to be 2x squared, minus b squared, minus 1 squared. 2x squared is 4x squared. 1 squared is just 1, so minus 1. So it's going to be 4x squared minus 1. Let's do one more of these, just to really hit the point home. I'll just focus on multiplication right now. If I have 5a minus 2b, and I'm multiplying that times 5a plus 2b. And remember, this only applies when I have at a product of a sum and a difference. That's the only time that I can use this. And I've shown you why. And if you're ever in doubt, just multiply it out. It'll take you a little bit longer. And you'll see the terms canceling out. You can't do this for just any binomial multiplication. You saw that earlier in the video, when we were multiplying, when we were taking squares. So this is going to be, using the pattern, it's going to be 5a squared minus 2b squared, which is equal to 25 a squared minus 4b squared. And, well, I'll leave it there, and I'll see you in the next video.


text:
## How to calculate classification accuracy with confusion matrix?

2

I have Train and Test data, how to calculate classification accuracy with confusion matrix ? Thanks

@attribute outlook {sunny, overcast, rainy}
@attribute temperature {hot, mild, cool}
@attribute humidity {high, normal}
@attribute windy {TRUE, FALSE}
@attribute play {yes, no}

Train:

1   sunny       hot     high    FALSE   no
2   sunny       hot     high    TRUE    no
3   overcast    hot     high    FALSE   yes
4   rainy       mild    high    FALSE   yes
5   rainy       cool    normal  FALSE   yes
6   rainy       cool    normal  TRUE    no
7   sunny       cool    normal  FALSE   yes
8   rainy       mild    normal  FALSE   yes
9   sunny       mild    normal  TRUE    yes
10  overcast    mild    high    TRUE    yes
11  overcast    hot     normal  FALSE   yes
12  rainy       mild    high    TRUE    no

Test:

overcast    cool    normal  TRUE    yes
sunny       mild    high    FALSE   no

Rules found:

(humidity,normal), (windy,FALSE) -> (play,yes) [Support=0.33 , Confidence=1.00 , Correctly Classify= 4, 8, 9, 12]
(outlook,overcast) -> (play,yes) [Support=0.25 , Confidence=1.00 , Correctly Classify= 2, 11]
(outlook,rainy), (windy,FALSE) -> (play,yes) [Support=0.25 , Confidence=1.00 , Correctly Classify= 3]
(outlook,sunny), (temperature,hot) -> (play,no) [Support=0.17 , Confidence=1.00 , Correctly Classify= 0, 1]
(outlook,sunny), (humidity,normal) -> (play,yes) [Support=0.17 , Confidence=1.00 , Correctly Classify= 10]
(outlook,rainy), (windy,TRUE) -> (play,no) [Support=0.17 , Confidence=1.00 , Correctly Classify= 5, 13]

How to calculate classification accuracy with confusion matrix ? Thanks – Xuan Dung – 2014-10-27T11:14:43.897

2

A confusion matrix is a cross tabulation of your predicted values against the true observed values, and (test) accuracy is emperical rate of correct predictions. So in this case you'll need to

1. Predict the 'play' attribute for your test set. (Currently you don't have a method to predict your second test case, so for the sake of argument let's assume your model would predict yes for the sunny example.
2. The following method of keeping track of your predictions is reffered to as a confusion matrix. The top labels are prediced
               Predicted
+----------------+
¦     ¦ yes ¦ no ¦
Oserved  ¦ yes ¦ 1   ¦ 1  ¦
¦ no  ¦ 0   ¦ 0  ¦
+----------------+

Here the first 1 is from your first test case and the second 1 is from the misclassified second test case.

1. Calculate accuracy,

Accuracy = (# correct predictions)/(# total predictions) = 1 / 2 = .50.

0

It is classify test objects: "In classification, let R be the set of generated rules and T the training data. The basic idea of the proposed method is to choose a set of high confidence rules in R to cover T. In classifying a test object, the first rule in the set of rules that matches the test object condition classifies it. This process ensures that only the highest ranked rules classify test objects."

Suppose 1 test case is (overcast, cool, normal, TRUE). Look through the rules top to bottom and see if any of the conditions are matched. The first rule for example tests the outlook feature. The value doesn't match, so the rule isn't matched. Move on to the next rule. And so on. In this case, rule 2 matches the test case and the classification for the play variable is "yes". The second test case is misclassified.

Thanks


text:
# Subjective video quality

Subjective video quality is video quality as experienced by humans. It is concerned with how video is perceived by a viewer (also called "observer" or "subject") and designates their opinion on a particular video sequence. It is related to the field of Quality of Experience. Measuring subjective video quality is necessary because objective quality assessment algorithms such as PSNR have been shown to correlate poorly with subjective ratings. Subjective ratings may also be used as ground truth to develop new algorithms.

Subjective video quality tests are psychophysical experiments in which a number of viewers rate a given set of stimuli. These tests are quite expensive in terms of time (preparation and running) and human resources and must therefore be carefully designed.

In subjective video quality tests, typically, SRCs ("Sources", i.e. original video sequences) are treated with various conditions (HRCs for "Hypothetical Reference Circuits") to generate PVSs ("Processed Video Sequences").[1]

## Measurement

The main idea of measuring subjective video quality is similar to the mean opinion score (MOS) evaluation for audio. To evaluate the subjective video quality of a video processing system, the following steps are typically taken:

• Choose original, unimpaired video sequences for testing
• Choose settings of the system that should be evaluated
• Apply settings to the SRC, which results in the test sequences
• Choose a test method, describing how sequences are presented to viewers and how their opinion is collected
• Invite a panel of viewers
• Carry out testing in a specific environment (e.g. a laboratory context) and present each PVS in a certain order to every viewer
• Calculate rating results for individual PVSs, SRCs and HRCs, e.g. the MOS

Many parameters of the viewing conditions may influence the results, such as room illumination, display type, brightness, contrast, resolution, viewing distance, and the age and educational level of viewers. It is therefore advised to report this information along with the obtained ratings.

### Source selection

Typically, a system should be tested with a representative number of different contents and content characteristics. For example, one may select excerpts from contents of different genres, such as action movies, news shows, and cartoons. The length of the source video depends on the purpose of the test, but typically, sequences of no less than 10 seconds are used.

The amount of motion and spatial detail should also cover a broad range. This ensures that the test contains sequences which are of different complexity.

Sources should be of pristine quality. There should be no visible coding artifacts or other properties that would lower the quality of the original sequence.

### Settings

The design of the HRCs depends on the system under study. Typically, multiple independent variables are introduced at this stage, and they are varied with a number of levels. For example, to test the quality of a video codec, independent variables may be the video encoding software, a target bitrate, and the target resolution of the processed sequence.

It is advised to select settings that result in ratings which cover the full quality range. In other words, assuming an Absolute Category Rating scale, the test should show sequences that viewers would rate from bad to excellent.

### Viewers

#### Number of viewers

Viewers are also called "observers" or "subjects". A certain minimum number of viewers should be invited to a study, since a larger number of subjects increases the reliability of the experiment outcome, for example by reducing the standard deviation of averaged ratings. Furthermore, there is a risk of having to exclude subjects for unreliable behavior during rating.

The minimum number of subjects that are required for a subjective video quality study is not strictly defined. According to ITU-T, any number between 4 and 40 is possible, where 4 is the absolute minimum for statistical reasons, and inviting more than 40 subjects has no added value. In general, at least 15 observers should participate in the experiment. They should not be directly involved in picture quality evaluation as part of their work and should not be experienced assessors.[2] In other documents, it is also claimed that at minimum 10 subjects are needed to obtain meaningful averaged ratings.[3]

However, most recommendations for the number of subjects have been designed for measuring video quality encountered by a home television or PC user, where the range and diversity of distortions tend to be limited (e.g., to encoding artifacts only). Given the large ranges and diversity of impairments that may occur on videos captured with mobile devices and/or transmitted over wireless networks, generally, a larger number of human subjects may be required.

Brunnström and Barkowsky have provided calculations for estimating the minimum number of subjects necessary based on existing subjective tests.[4] They claim that in order to ensure statistically significant differences when comparing ratings, a larger number of subjects than usually recommended may be needed.

#### Viewer selection

Viewers should be non-experts in the sense of not being professionals in the field of video coding or related domains. This requirement is introduced to avoid potential subject bias.[2]

Typically, viewers are screened for normal vision or corrected-to-normal vision using Snellen charts. Color blindness is often tested with Ishihara plates.[2]

There is an ongoing discussion in the QoE community as to whether a viewer's cultural, social, or economic background has a significant impact on the obtained subjective video quality results. A systematic study involving six laboratories in four countries found no statistically significant impact of subject's language and culture / country of origin on video quality ratings.[5]

### Test environment

Subjective quality tests can be done in any environment. However, due to possible influence factors from heterogenous contexts, it is typically advised to perform tests in a neutral environment, such as a dedicated laboratory room. Such a room may be sound-proofed, with walls painted in neutral grey, and using properly calibrated light sources. Several recommendations specify these conditions.[6][7] Controlled environments have been shown to result in lower variability in the obtained scores.[5]

#### Crowdsourcing

Crowdsourcing has recently been used for subjective video quality evaluation, and more generally, in the context of Quality of Experience.[8] Here, viewers give ratings using their own computer, at home, rather than taking part in a subjective quality test in laboratory rooms. While this method allows for obtaining more results than in traditional subjective tests at lower costs, the validity and reliability of the gathered responses must be carefully checked.[9]

The world's largest crowdsourced picture quality database was made publicly available as the LIVE In the Wild Image Quality Challenge Database. It contains more than 350,000 human quality judgments from more than 8,000 human subjects. The human judgments were made on more than 1,100 mobile camera photographs afflicted by varying degrees of authentic distortions, instead of synthetically-introduced distortions as has been done in nearly every prior database.

### Analysis of results

Opinions of viewers are typically averaged into the mean opinion score (MOS). To this aim, the labels of categorical scales may be translated into numbers. For example, the responses "bad" to "excellent" can be mapped to the values 1 to 5, and then averaged. MOS values should always be reported with their statistical confidence intervals so that the general agreement between observers can be evaluated.

#### Subject screening

Often, additional measures are taken before evaluating the results. Subject screening is a process in which viewers whose ratings are considered invalid or unreliable are rejected from further analysis. Invalid ratings are hard to detect, as subjects may have rated without looking at a video, or cheat during the test. The overall reliability of a subject can be determined by various procedures, some of which are outlined in ITU-R and ITU-T recommendations.[2][7] For example, the correlation between a person's individual scores and the overall MOS, evaluated for all sequences, is a good indicator of their reliability in comparison with the remaining test participants.

While rating stimuli, humans are subject to biases. These may lead to different and inaccurate scoring behavior and consequently result in MOS values that are not representative of the “true quality” of a stimulus. In the recent years, advanced models have been proposed that aim at formally describing the rating process and subsequently recovering noisiness in subjective ratings. According to Janowski et al., subjects may have an opinion bias that generally shifts their scores, as well as a scoring imprecision that is dependent on the subject and stimulus to be rated.[10] Li et al. have proposed to differentiate between subject inconsistency and content ambiguity.[11]

## Standardized testing methods

There are many ways to select proper sequences, system settings, and test methodologies. A few of them have been standardized. They are thoroughly described in several ITU-R and ITU-T recommendations, among those ITU-R BT.500[7] and ITU-T P.910.[2] While there is an overlap in certain aspects, the BT.500 recommendation has its roots in broadcasting, whereas P.910 focuses on multimedia content.

A standardized testing method usually describes the following aspects:

• how long an experiment session lasts
• where the experiment takes place
• how many times and in which order each PVS should be viewed
• whether ratings are taken once per stimulus (e.g. after presentation) or continuously
• whether ratings are absolute, i.e. referring to one stimulus only, or relative (comparing two or more stimuli)
• which scale ratings are taken on

Another recommendation, ITU-T P.913,[6] gives researchers more freedom to conduct subjective quality tests in environments different from a typical testing laboratory, while still requiring them to report all details necessary to make such tests reproducible.

### Examples

Below, some examples of standardized testing procedures are explained.

#### Single-Stimulus

• ACR (Absolute Category Rating):[2] each sequence is rated individually on the ACR scale. The labels on the scale are "bad", "poor", "fair", "good", and "excellent", and they are translated to the values 1, 2, 3, 4 and 5 when calculating the MOS.
• ACR-HR (Absolute Category Rating with Hidden Reference): a variation of ACR, in which an original unimpaired source sequence is shown in addition to the impaired sequences, without informing the subjects of its presence (hence, "hidden"). The ratings are calculated as differential scores between the reference and the impaired versions. The differential score is defined as the score of the PVS minus the score given to the hidden reference, plus the number of points on the scale. For example, if a PVS is rated as “poor", and its corresponding hidden reference as “good", then the rating is ${\textstyle 2-4+5=3}$ . When these ratings are averaged, the result is not a MOS, but a differential MOS ("DMOS").
• SSCQE (Single Stimulus Continuous Quality Rating):[7] a longer sequence is rated continuously over time using a slider device (a variation of a fader), on which subjects rate the current quality. Samples are taken in regular intervals, resulting in a quality curve over time rather than a single quality rating.

#### Double-stimulus or multiple stimulus

• DSCQS (Double Stimulus Continuous Quality Scale):[7] the viewer sees an unimpaired reference and the impaired sequence in a random order. They are allowed to re-view the sequences, and then rate the quality for both on a continuous scale labeled with the ACR categories.
• DSIS (Double Stimulus Impairment Scale)[7] and DCR (Degradation Category Rating):[2] both refer to the same method. The viewer sees an unimpaired reference video, then the same video impaired, and after that they are asked to vote on the second video using a so-called impairment scale (from "impairments are imperceptible" to "impairments are very annoying").
• PC (Pair Comparison):[2] instead of comparing an unimpaired and impaired sequence, different impairment types (HRCs) are compared. All possible combinations of HRCs should be evaluated.

### Choice of methodology

Which method to choose largely depends on the purpose of the test and possible constraints in time and other resources. Some methods may have fewer context effects (i.e. where the order of stimuli influences the results), which are unwanted test biases.[12] In ITU-T P.910, it is noted that methods such as DCR should be used for testing the fidelity of transmission, especially in high quality systems. ACR and ACR-HR are better suited for qualification tests and – due to giving absolute results – comparison of systems. The PC method has a high discriminatory power, but it requires longer test sessions.

## Databases

The results of subjective quality tests, including the used stimuli, are called databases. A number of subjective picture and video quality databases based on such studies have been made publicly available by research institutes. These databases – some of which have become de facto standards – are used globally by television, cinematic, and video engineers around the world to design and test objective quality models, since the developed models can be trained against the obtained subjective data.

Examples of these databases include the popular collection of public-domain subjective picture quality databases available at the Laboratory for Image & Video Engineering (LIVE) website. The LIVE databases were the first high quality, large subjective picture quality databases to be made available to the public free of charge and without withholding of any data portion. Recent work includes the databases by the Institut de Recherche en Communications et Cybernétique de Nantes (IRCCyN). A list of European databases related to the field of Quality of Experience can be found in the QUALINET databases.

## References

1. ^
2. ^ Winkler, Stefan. "On the properties of subjectiveratings in video quality experiments". Proc. Quality of Multimedia Experience, 2009.
3. ^ Brunnström, Kjell; Barkowsky, Marcus (2018-09-25). "Statistical quality of experience analysis: on planning the sample size and statistical significance testing". Journal of Electronic Imaging. 27 (5): 053013. Bibcode:2018JEI....27e3013B. doi:10.1117/1.jei.27.5.053013. ISSN 1017-9909. S2CID 53058660.
4. ^ a b Pinson, M. H.; Janowski, L.; Pepion, R.; Huynh-Thu, Q.; Schmidmer, C.; Corriveau, P.; Younkin, A.; Callet, P. Le; Barkowsky, M. (October 2012). "The Influence of Subjects and Environment on Audiovisual Subjective Tests: An International Study" (PDF). IEEE Journal of Selected Topics in Signal Processing. 6 (6): 640–651. Bibcode:2012ISTSP...6..640P. doi:10.1109/jstsp.2012.2215306. ISSN 1932-4553. S2CID 10667847.
5. ^ a b
6. ^ Hossfeld, Tobias (2014-01-15). "Best Practices for QoE Crowdtesting: QoE Assessment With Crowdsourcing". IEEE Transactions on Multimedia. 16 (2): 541–558. doi:10.1109/TMM.2013.2291663. S2CID 16862362.
7. ^ Hossfeld, Tobias; Hirth, Matthias; Redi, Judith; Mazza, Filippo; Korshunov, Pavel; Naderi, Babak; Seufert, Michael; Gardlo, Bruno; Egger, Sebastian (October 2014). "Best Practices and Recommendations for Crowdsourced QoE - Lessons learned from the Qualinet Task Force "Crowdsourcing"". hal-01078761. Cite journal requires `|journal=` (help)
8. ^ Janowski, Lucjan; Pinson, Margaret (2015). "The Accuracy of Subjects in a Quality Experiment: A Theoretical Subject Model". IEEE Transactions on Multimedia. 17 (12): 2210–2224. doi:10.1109/tmm.2015.2484963. ISSN 1520-9210. S2CID 22343847.
9. ^ Li, Zhi; Bampis, Christos G. (2017). "Recover Subjective Quality Scores from Noisy Measurements". 2017 Data Compression Conference (DCC). IEEE: 52–61. arXiv:1611.01715. doi:10.1109/dcc.2017.26. ISBN 9781509067213. S2CID 14251604.
10. ^ Pinson, Margaret and Wolf, Stephen. "Comparing Subjective Video Quality Testing Methodologies". SPIE Video Communications and Image Processing Conference, Lugano, Switzerland, July 2003.