The 1D wave equation is $$\frac{\partial^2 u(x,t)}{\partial t^2} = c^2 \frac{\partial^2 u(x,t)}{\partial x^2}$$ where $c$ is the wave speed, $c^2=E/\rho$, $E$ is the Young's modulus and $\rho$ is the density of the material. See the link. I want to analyze the physics in this video of longitudinal wave. When a mass impacts one end of a long rod with the material properties of $E$ and $\rho$ at $x=0$, longitudinal waves will be generated. Please note the $\textbf{initial conditions}$: $$v(x,t)=\dot{u}(x,t)=0$$ $$u(x,0)=0$$ where $v(x,t)$ represents the particle velocity in the medium. The $u(x,t)$ should be zero at time $t=0$.

Note that for any given time $0<t_1< L/c$, the longitudinal waves will travel to $$x_1=c\cdot t_1$$ For $x\geq x_1$, it is required that $$v(x,t_1)=0$$ $$u(x,t_1)=0$$ The waves will travel through the rod to the other end. Assume the length of the rod is $L$. The wave speed is $c$. At time $t=L/c$ the waves will reach $x=L$. Depending on the conditions and the other end, the waves can change and then be reflected.

  1. Free end condition. That is to say, there is nothing in contact with the other end (That is $x=L$). The particle in the medium can move in the free end. What will happen when the wave reaches there?
  2. Fixed end condition. In this condition, the particles in the medium at the fixed end cannot move. The particle velocity in the fixed end is strictly zero when the waves travel $x=L$. What will happen in the fixed end?

How should I analyze the two case scenarios for the 1D wave equation in Mathematica? I want to see how $u(x,t)$ will vary over $t$ and $x$.

Thank you very much!

Edited on 10/07/2022: Using the solution from Ulrich Neumann, the following animation presents the results.I have a question. Why will the $u(x,t)$ increase over time? I did not expect that.

Solution by Ulrich

  • $\begingroup$ I'm afraid your b.c. at $x=0$ isn't a proper modeling for the scene in the video, see this post for more info: mathematica.stackexchange.com/a/71945/1871 Something like $u(0,t)=\sin ^{30}(t)$ (this is just an example) may be a better choice. $\endgroup$
    – xzczd
    Oct 6, 2022 at 3:59
  • $\begingroup$ Have a look here for various boundary conditions and examples for the wave equation $\endgroup$
    – user21
    Oct 6, 2022 at 4:01
  • $\begingroup$ You significantly changed your question without notification! My answer still describes the wave for the original problem u[0,t]==0, Derivative[0,1][u][x,t]==0 and impact at x==0! $\endgroup$ Oct 7, 2022 at 6:44
  • $\begingroup$ @UlrichNeumann Thanks for your help. Is you current solution applicable for the free end condition? What will be appropriate boundary condition for fixed end condition? $\endgroup$
    – fhk
    Oct 7, 2022 at 14:29
  • $\begingroup$ @fhk Yes it is with free end condition. Inside NDSolveValu[...] you'll find the condition (*u[t,1]\[Equal]0,*) for fixed end. Remove comment and you get solution for fixed end. $\endgroup$ Oct 7, 2022 at 14:42

2 Answers 2


If I understand your question right, you are considering a rod, with an

impact Derivative[0,1][u][t,0]==DiracDelte[t] at x==0


boundary condition u[t,1]==0 or Derivative[0,1][u][t,1]==0 at x==1!

Analytical solution might be described using GreenFunction , unfortunately Mathematica doesn't solve.

Remains the numerical solution:

The impact should be modelled by a numerical approximation UnitTriangle[t/(eps/2) - 1] (2/eps),eps<<1 of the DiracDelta -function. Additionally I introduced a small damping to smooth simulation results.

For the free end condition Derivative[ 0, 1][u][t, 1]==0 try

eps = 1/10 ;
U = NDSolveValue[{Derivative[2, 0][u][t, x] -Derivative[ 0, 2][u][t, x] + 2 /100 Derivative[1, 0][u][t, x] == 
UnitTriangle[t/(eps/2) - 1] 2/eps NeumannValue[1, x ==0],(*u[t,1]\[Equal]0,*)Derivative[1, 0][u][0, x] == 0, u[0, x] == 0}, 
u, {t, 0, 5}, {x, 0, 1}   , 
Method -> {"MethodOfLines", "TemporalVariable" -> t,"SpatialDiscretization" -> {"FiniteElement" ,"MeshOptions"-> { "MaxCellMeasure" -> 1/50} }}  , MaxStepSize -> eps/10]

Plot3D[Evaluate[ U[t, x]], {x, 0, 1}, {t, 0, 5} , PlotPoints ->100,MeshFunctions -> {#2 &}, PlotRange -> All]

enter image description here

Plot shows the waves for differnt times and also the global movement of the free rod !

Add u[t,1]==0 for the second end condition.

addendum fixed end

eps = 1/10 ; tsim = 5;
Uf = NDSolveValue[{Derivative[2, 0][u][t, x] - 
     Derivative[ 0, 2][u][t, x] + 
     2  10^-2 Derivative[1, 0][u][t, x] == 
    UnitTriangle[t/(eps/2) - 1] 2/eps NeumannValue[1, x == 0],  
   u[t, 1] == 0,   
    Derivative[1, 0][u][0, x] == 0, u[0, x] == 0}, 
  u, {t, 0, tsim}, {x, 0, 1}   , 
  Method -> {"MethodOfLines", "TemporalVariable" -> t, 
    "SpatialDiscretization" -> {"FiniteElement" , 
      "MeshOptions" -> { "MaxCellMeasure" -> 1/50} }}  , 
  MaxStepSize -> eps/10]

  list = Table[Plot[Uf[t, x], {x, 0, 1}, PlotRange -> {-1, 1}], {t,Subdivide[0, tsim, 100]}];

enter image description here

  • $\begingroup$ Excellent! Can you explain the coefficients in the PDE? In NDSolveValue[], you define the PDE. 1. What does "2/100" mean in the PDE? 2. If I want to specify an initial particle velocity at x==0, how should I do? 3. Why does we still need some coefficient in the right hand side of the PDE? $\endgroup$
    – fhk
    Oct 7, 2022 at 16:46
  • $\begingroup$ @fhk Righthandside describes the impact and 2 /100 Derivative[1, 0][u][t, x] is a small damping $\endgroup$ Oct 7, 2022 at 18:36
  • $\begingroup$ @fhk You need to know velocity and mass of the particle $\endgroup$ Oct 7, 2022 at 18:37
  • $\begingroup$ The boundary condition for free end $u_x(t,x)=0$ was not enforced at $x=1$. I do not see '''Derivative[0,1][u][t,1]==0 at x==1!''' in your code. It still produces the results. $\endgroup$
    – fhk
    Oct 8, 2022 at 2:39
  • $\begingroup$ @fhk Condition for free end is implicitly included in the first part of my answer. You could add an additional NeumannValue[0, x ==1] $\endgroup$ Oct 8, 2022 at 9:06

This is an incorrect solution for the fixed end condition where $u(L,t)=0$:

define some parameters:

c = 500;
L = 5;
t0 = 0.1;
U = NDSolveValue[{
   D[u[t, x], {t, 2}] == c^2 D[u[t, x], {x, 2}],
   u[0, x] == 0,
   u[t, 0] == Sin[80 \[Pi] t],
   u[t, L] == 0},
  u, {t, 0, t0}, {x, 0, L}]

Show the results:

list = Table[
   Plot[U[t, x], {x, 0, L}, PlotRange -> {-1, 1}], {t, 0, t0, 0.003}];

enter image description here

  • $\begingroup$ An impact at x==0&&t==0 causes a wave, starting from x==0 and travelling in positive x-direction. Your "solution" doesn't show this effect! $\endgroup$ Oct 7, 2022 at 7:11
  • $\begingroup$ I agree well. Please propose some initial conditions or boundary conditions. How should we fix this problem? $\endgroup$
    – fhk
    Oct 7, 2022 at 13:45
  • $\begingroup$ See my answer including the Boundaryconditions and NeumannValue ! $\endgroup$ Oct 7, 2022 at 14:13
  • $\begingroup$ As it’s currently written, your answer is unclear. Please edit to add additional details that will help others understand how this addresses the question asked. You can find more information on how to write good answers in the help center. $\endgroup$
    – Community Bot
    Oct 8, 2022 at 2:56

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.