Как построить экспериментальные данные в особых точках?

Question

У меня есть набор экспериментальных данных (P), и я хочу получить график "экспериментальное против предсказанного".Чтобы сделать это, я использую другой набор данных, которые зависят от P (Q), строю график ScatterPlot, использую подходящее соответствие, затем получаю линию регрессии и использую ее коэффициенты в соответствующем дифференциальном уравнении.График P выглядит хорошо, но мне нужно добавить туда экспериментальные данные.Для простоты я использовал интервал t=0..150.

Как я могу построить экспериментальные данные так, чтобы P(0) = Pvals[1], P(10)=Pvals[2] и т. Д.?Кроме того, как я могу распространять данные (скажем, у меня есть t=0..800 и я хочу построить Pvals так, чтобы P(0) = Pvals[1] and P(800) = Pvals[16])?

Pvals := [3.929, 5.308, 7.24, 9.638, 12.866, 17.069, 23.192, 31.433, 38.558, 50.156, 62.948, 
75.996, 91.972, 105.711, 122.775, 131.669]:

for i to 15 do Qval[i] := .1*(Pvals[i+1]/Pvals[i]-1); end do:
Qvals := [seq(Qval[i], i = 1 .. 15), 0.144513895e-1]:               
    with(Statistics);
ScatterPlot(Pvals, Qvals, fit = [a*v^2+b*v+c, v], thickness = 3, 
legend = [points = "Point data", fit = typeset("fit to a", 2^nd, "degree polynomial")]);

with(CurveFitting);
LeastSquares(Pvals, Qvals, v, curve = a*v^2+b*v+c);

de := diff(P(t), t) = (0.370152282598477e-1-0.272504103112702e-3*P(t))*P(t);

sol := dsolve({de, P(0) = 3.929}, P(t));

P := plot(rhs(sol), t = 0 .. 160);

Answer 1

Я не уверен, что полностью следую вашей методологии.Но это что-то вроде того, что вы пытаетесь выполнить?

restart;
with(Statistics):

Pvals := [3.929, 5.308, 7.24, 9.638, 12.866, 17.069, 23.192, 31.433,
          38.558, 50.156, 62.948, 75.996, 91.972, 105.711, 122.775, 131.669]:

for i to 15 do Qval[i] := .1*(Pvals[i+1]/Pvals[i]-1); end do:
Qvals := [seq(Qval[i], i = 1 .. 15), 0.144513895e-1]:

form := a*v^2+b*v+c:

CF := CurveFitting:-LeastSquares(Pvals, Qvals, v, curve = form);

       CF := 0.0370152282598477 - 0.000272504103112702 v

                               -7  2
              + 5.60958249026713 10   v

Теперь я использую CF в DE (поскольку я не понимаю, почему вы отбросили термин v ^ 2),

#de := diff(P(t), t) = (0.370152282598477e-1-0.272504103112702e-3*P(t))*P(t);
de := diff(P(t), t) = eval(CF, v=P(t))*P(t);

           d         /                                              
    de := --- P(t) = \0.0370152282598477 - 0.000272504103112702 P(t)
           dt                                                       

                            -7     2\     
       + 5.60958249026713 10   P(t) / P(t)

Я буду использовать опцию numeric команды dsolve и получу процедуру, которая вычисляет P(t) для числовых t значений.

sol := dsolve({de, P(0) = 3.929}, P(t), numeric, output=listprocedure ):

Pfunc := eval(P(t), sol);

              Pfunc := proc(t)  ...  end;

Pfunc(0.0), Pvals[1];

                3.92900000000000, 3.929

Теперь немного изменим масштаб (чтоопять же, мое предположение о вашей цели),

endpt := fsolve(Pfunc(t)-Pvals[16]);

                  endpt := 135.2246055

Pfunc(endpt), Pvals[16];

               131.669000003321, 131.669

plot(Pfunc(t), t=0 .. endpt, size=[500,200]);

a,b,N := 0.0, 800.0, nops(Pvals);

                a, b, N := 0., 800.0, 16

Pfuncscaled := proc(t) 
                 if not t::numeric then
                   return 'procname'(args);
                 end if;
                 Pfunc(t*endpt/b);
               end proc:

Pfuncscaled(0), Pvals[1];

                3.92900000000000, 3.929

Pfuncscaled(800), Pvals[N];

               131.669000003321, 131.669

PLscaled := plot( Pfuncscaled(t), t=a .. b,
                  color=red, size=[500,200] );

Теперь, чтобы отобразить Pdata и 0 .. 800,

V := Vector(N, (i)->a+(i-1)*(b-a)/(N-1)):

V[1], V[-1];

                    0., 800.0000000

Pdatascaled := plot( < V | Vector(Pvals) >,
                     color=blue, size=[500,200],
                     style=pointline, symbol=solidcircle );

И, отображая измененные данные вместе спроцедура изменения масштаба от dsolve,

plots:-display( PLscaled, Pdatascaled, size=[500,500] );