Solving a damped wave equation

Question

I am trying to solve the equation

$$ \frac{d^2u}{dt^2}-\frac{d^2}{dx^2}\left(c_s^2u+\nu\frac{du}{dt}\right)=0 $$

with initial conditions

$$u(x, 0)=0$$

$$\frac{du}{dt}|_{t=0}=0$$

and boundary conditions

$$ {\frac{du}{dx}}|_{x=0,1}=a\sin{(\omega_d t)}-b\cos{(\omega_d t)}$$

My attempts so far are

ClearAll[u, x, t, a, b, c, w, n];

c = 1;
n = 1;
a = 1;
b = 1;
w = 1;

pde = D[u[t, x], t, t] - c*D[u[t, x], x, x] - n*D[u[t, x], x, x, t] == 0;
ics = {u[0,x]==0};
bcs = 
  {(D[u[t,x],x] /. x->0) == a*Sin[w*t]-b*Cos[w*t], 
   (D[u[t,x],x] /. x->1) == a*Sin[w*t]-b*Cos[w*t]};

sol = NDSolveValue[{pde, ics, bcs}, u, {x, 0, 1},{t, 0, 10}]

but I am receiving multiple error messages. Among them:

NDSolveValue::fembdnl: The dependent variable in (u^(0,1))[t,0]==-Cos[t]+Sin[t] in the boundary condition DirichletCondition[(u^(0,1))[t,0]==-Cos[t]+Sin[t],x==0.] needs to be linear.

NDSolveValue::fembdnl: The dependent variable in (u^(0,1))[t,0]==-Cos[t]+Sin[t] in the boundary condition DirichletCondition[(u^(0,1))[t,0]==-Cos[t]+Sin[t],x==0.] needs to be linear.

NDSolveValue::femcmsd: The spatial derivative order of the PDE may not exceed two.

How can I solve it?

Answer 1

It is necessary to divide the equation into two (second order), add the initial data and boundary conditions. Here is the code without adding boundary conditions

ClearAll[u, x, t, a, b, c, w, n];

c = 1;
n = 1;
a = 1;
b = 1;
w = 1;

pde = {D[u[t, x], t, t] - c*D[u[t, x], x, x] - n*D[v[t, x], t] == 0, 
   v[t, x] == D[u[t, x], x, x]};
ics = {u[0, x] == 0, v[0, x] == 0, Derivative[1, 0][u][0, x] == 0};
bcs = {(D[u[t, x], x] /. x -> 0) == 
    If[t <= 10^-6, 0, 
     a*Sin[w*t] - b*Cos[w*t]], (D[u[t, x], x] /. x -> 1) == 
    If[t <= 10^-6, 0, a*Sin[w*t] - b*Cos[w*t]]};

{U, V} = NDSolveValue[{pde, ics, bcs}, {u, v}, {x, 0, 1}, {t, 0, 10}, 
  Method -> {"IndexReduction" -> Automatic, 
    "EquationSimplification" -> "Residual", 
    "PDEDiscretization" -> {"MethodOfLines", 
      "SpatialDiscretization" -> {"TensorProductGrid", 
        "MinPoints" -> 141, "MaxPoints" -> 141, 
        "DifferenceOrder" -> 2}}}]

DensityPlot[U[t, x], {x, 0, 1}, {t, 0, 10}, 
 ColorFunction -> "Rainbow", PlotLegends -> Automatic, 
 FrameLabel -> Automatic]

We can use boundary conditions twice, then there are no warnings

ClearAll[u, x, t, a, b, c, w, n];

c = 1;
n = 1;
a = 1;
b = 1;
w = 1;
f[t_] := If[t < 10^-6, 0, 1];
pde = {D[v[t, x], t] - c*D[u[t, x], x, x] - n*D[v[t, x], x, x] == 0, 
   v[t, x] == D[u[t, x], t]};
ics = {u[0, x] == 0, v[0, x] == 0};
bcs = {(D[u[t, x], x] /. x -> 0) == (a*Sin[w*t] - b*Cos[w*t]) f[
      t], (D[u[t, x], x] /. x -> 1) == a*Sin[w*t] - b*Cos[w*t] f[t]};
bcs1 = {(D[v[t, x], x] /. x -> 0) == 
    w (a*Cos[w*t] + b*Sin[w*t]) f[t], (D[v[t, x], x] /. x -> 1) == 
    w (a*Cos[w*t] + b*Sin[w*t]) f[t]};

{U, V} = NDSolveValue[{pde, ics, bcs, bcs1}, {u, v}, {x, 0, 1}, {t, 0,
     10}, Method -> {"IndexReduction" -> Automatic, 
     "EquationSimplification" -> "Residual", 
     "PDEDiscretization" -> {"MethodOfLines", 
       "SpatialDiscretization" -> {"TensorProductGrid", 
         "MinPoints" -> 141, "MaxPoints" -> 141, 
         "DifferenceOrder" -> 2}}}];

{DensityPlot[U[t, x], {x, 0, 1}, {t, 0, 10}, 
  ColorFunction -> "Rainbow", PlotLegends -> Automatic, 
  FrameLabel -> Automatic, PlotLabel -> "u"], 
 DensityPlot[V[t, x], {x, 0, 1}, {t, 0, 10}, 
  ColorFunction -> "Rainbow", PlotLegends -> Automatic, 
  FrameLabel -> Automatic, PlotLabel -> "v"]}