我要求解一维波动方程，一维杆左端（x=0）固定，右端（x=L）自由，初始时刻杆处于静力平衡状态，在右端（x=L）加载一个自己定义的荷载P(t)。我的代码如下
EE = 30000000000;
\[Rho] = 2600;
\[Mu] = 0.2;
a = Sqrt[EE/\[Rho]*(1 - \[Mu])/((1 - 2*\[Mu])*(1 + \[Mu]))];
\[Xi] = 1500;
Pmax = 500000000;
Pmid = 500000000/2;
t1 = 0.01/1000;
t2 = 0.1/1000;
t3 = 5/1000;
Pt = \!\(\*
TagBox[GridBox[{
{"\[Piecewise]", GridBox[{
{
RowBox[{
RowBox[{"Pmax", "/", "t1"}], "*", "t"}],
RowBox[{"0", "<=", "t", "<=", "t1"}]},
{
RowBox[{"Pmax", "-",
RowBox[{
RowBox[{
RowBox[{"(",
RowBox[{"Pmax", "-", "Pmid"}], ")"}], "/",
RowBox[{"(",
RowBox[{"t2", "-", "t1"}], ")"}]}], "*",
RowBox[{"(",
RowBox[{"t", "-", "t1"}], ")"}]}]}],
RowBox[{"t1", "<=", "t", "<=", "t2"}]},
{
RowBox[{"Pmid", "-",
RowBox[{
RowBox[{"Pmid", "/",
RowBox[{"(",
RowBox[{"t3", "-", "t2"}], ")"}]}], "*",
RowBox[{"(",
RowBox[{"t", "-", "t2"}], ")"}]}]}],
RowBox[{"t2", "<=", "t", "<=", "t3"}]},
{"0",
RowBox[{"t", ">=", "t3"}]}
},
AllowedDimensions->{2, Automatic},
Editable->True,
GridBoxAlignment->{"Columns" -> {{Left}}, "Rows" -> {{Baseline}}},
GridBoxItemSize->{"Columns" -> {{Automatic}}, "Rows" -> {{1.}}},
GridBoxSpacings->{"Columns" -> {
Offset[0.27999999999999997`], {
Offset[0.84]},
Offset[0.27999999999999997`]}, "Rows" -> {
Offset[0.2], {
Offset[0.4]},
Offset[0.2]}},
Selectable->True]}
},
GridBoxAlignment->{"Columns" -> {{Left}}, "Rows" -> {{Baseline}}},
GridBoxItemSize->{"Columns" -> {{Automatic}}, "Rows" -> {{1.}}},
GridBoxSpacings->{"Columns" -> {
Offset[0.27999999999999997`], {
Offset[0.35]},
Offset[0.27999999999999997`]}, "Rows" -> {
Offset[0.2], {
Offset[0.4]},
Offset[0.2]}}],
"Piecewise",
DeleteWithContents->True,
Editable->False,
SelectWithContents->True,
Selectable->False]\);
L = 10;
s1 = NDSolve[{D[u[x, t], t, t] -
a^2 D[u[x, t], x, x] + \[Xi]*D[u[x, t], t] == 0, u[0, t] == 0,
Derivative[0, 1][u][x, 0] == 0, u[x, 0] == 0,
EE*Derivative[1, 0][u][L, t] == Pt}, u, {x, 0, L}, {t, 0, 20/1000}]
然后得到了两个提示
1.Warning: boundary and initial conditions are inconsistent.
2.Encountered non-numerical value for a derivative at t == 0.`.
请大家帮忙看一下

1. 边界条件和初始条件不匹配，你要改一下定解条件。
2. 分段函数的表达方式不是你以为的那样。用Piecewise写分段函数。
Clear["Global`*"]
EE = 30000000000;
\[Rho] = 2600;
\[Mu] = 0.2;
a = Sqrt[EE/\[Rho]*(1 - \[Mu])/((1 - 2*\[Mu])*(1 + \[Mu]))];
\[Xi] = 1500;
Pmax = 500000000;
Pmid = 500000000/2;
t1 = 0.01/1000;
t2 = 0.1/1000;
t3 = 5/1000;
Pt[t_] :=
Piecewise[{{Pmax/t1*t,
0 <= t <= t1}, {Pmax - (Pmax - Pmid)/(t2 - t1)*(t - t1),
t1 <= t <= t2}, {Pmid - Pmid/(t3 - t2)*(t - t2),
t2 <= t <= t3}, {0, t >= t3}}];
L = 10;
s1 = NDSolve[{D[u[x, t], t, t] -
a^2 D[u[x, t], x, x] + \[Xi]*D[u[x, t], t] == 0, u[0, t] == 0,
Derivative[0, 1][u][x, 0] == 0, u[x, 0] == 0,
EE*Derivative[1, 0][u][L, t] == Pt}, u, {x, 0, L}, {t, 0, 20/1000}]

感谢你的解答，我还有点疑问。我复制你的代码到新的nb文件后，在右端x=L的边界条件右端是Pt，可以计算出数值解u[x,t]，同时第2个提示消失了，只有第1个提示，但是用Plot3D画出来是一个值等于0的平面；我把Pt改为Pt[t]，重新计算，就出现我原问题中同样的1和2提示。
所以我想问，1.按你那样延时定义的方法，边界条件中应该是Pt还是P[t]呢？
2.关于第2个提示，“Encountered non-numerical value for a derivative at t == 0.`.”，这个应该是个错误吧，这是因为什么原因出现的呢？以及该怎么避免。
感谢你的帮助

……你撞上bug了。ibcinc警告很多时候不是问题，但是，如果定义问题时用了分段函数，它经常会触发ndnum警告，这是一个已经被Wolfram客服证实的bug，你不妨把你这个也报上去，以便引起他们的重视。更多内容可以参看Stackexchange的[discontinuity]标签下的问题。
对于你这个问题，可以通过使用连续函数近似原本的分段函数来规避问题。注意NDSolve的选项也要做相应调整：
mol[n:_Integer|{_Integer..}, o_:"Pseudospectral"] := {"MethodOfLines",
"SpatialDiscretization" -> {"TensorProductGrid", "MaxPoints" -> n,
"MinPoints" -> n, "DifferenceOrder" -> o}}
mol[tf:False|True,sf_:Automatic]:={"MethodOfLines",
"DifferentiateBoundaryConditions"->{tf,"ScaleFactor"->sf}}
appro = With[{k = 10^6}, ArcTan[k #]/Pi + 1/2 &];
EE = 30000000000;
\[Rho] = 2600;
\[Mu] = 0.2;
a = Sqrt[EE/\[Rho]*(1 - \[Mu])/((1 - 2*\[Mu])*(1 + \[Mu]))];
\[Xi] = 1500;
Pmax = 500000000;
Pmid = 500000000/2;
t1 = 0.01/1000;
t2 = 0.1/1000;
t3 = 5/1000;
Pt = Simplify`PWToUnitStep@ Piecewise[{{(Pmax/t1)*t, 0 <= t <= t1},
{Pmax - ((Pmax - Pmid)/(t2 - t1))*(t - t1), t1 <= t <= t2},
{Pmid - (Pmid/(t3 - t2))*(t - t2), t2 <= t <= t3}, {0, t >= t3}}] /.
UnitStep -> appro;
L = 10;
s1 = NDSolveValue[{D[u[x, t], t, t] -
a^2 D[u[x, t], x, x] + \[Xi]*D[u[x, t], t] == 0, u[0, t] == 0,
Derivative[0, 1][u][x, 0] == 0, u[x, 0] == 0,
EE*Derivative[1, 0][u][L, t] == Pt}, u, {x, 0, L}, {t, 0, 20/1000},
MaxSteps -> Infinity, Method -> Union[mol[25, 4], mol[True, 10^6]]]
检查一下边界：
Plot[D[s1[x, t], x] /. x -> L // Evaluate, {t, 0, 0.001}, PlotRange -> All]
嗯，近似度还是可以的。看样子取10^6是够的。
Plot3D[s1[x, t], {t, 0, 20/1000}, {x, 0, L}, PlotPoints -> 50]

指定纽曼条件和狄利克雷条件就可以了，但这种黑箱操作，很容易让人怀疑结果准不准，建议用多种方法验证一下。例如这里u的最大值为u=0.0294629，取在x=10,t=0.00292567处。
代码为，
Clear["Global`*"];
EE = 30000000000;
\[Rho] = 2600;
\[Mu] = 0.2;
a = Sqrt[EE/\[Rho]*(1 - \[Mu])/((1 - 2*\[Mu])*(1 + \[Mu]))];
\[Xi] = 1500;
Pmax = 500000000;
Pmid = 500000000/2;
t1 = 0.01/1000;
t2 = 0.1/1000;
t3 = 5/1000;
Pt[t_] = Piecewise[{{(Pmax/t1)*t, 0 <= t <= t1},
{Pmax - ((Pmax - Pmid)/(t2 - t1))*(t - t1), t1 <= t <= t2},
{Pmid - (Pmid/(t3 - t2))*(t - t2), t2 <= t <= t3}, {0, t >= t3}}];
L = 10;
sol = NDSolveValue[{D[u[x, t], t, t] -
a^2 D[u[x, t], x, x] + \[Xi]*D[u[x, t], t] ==
NeumannValue[Pt[t]/EE*a^2, x == L],
DirichletCondition[u[x, t] == 0, x == 0],
Derivative[0, 1][u][x, 0] == 0, u[x, 0] == 0}, u, {x, 0, L},
{t, 0, 20/1000}, Method -> "MethodOfLines"];
{Plot3D[sol[x, t], {t, 0, 20/1000}, {x, 0, L}, PlotPoints -> 50,
PlotTheme -> "Classic"],
Plot[sol[10, t], {t, 0, 20/1000}, Mesh -> Full, Frame -> True]}
Maximize[{sol[10, t], 0 <= t <= 20/1000}, t]