## 13821 Reputation

14 years, 211 days

## Not sure you want to "display" this,...

because it is a matrix with 702 entries. The attached will display the first 10 or so - al the others are accessible. If you really wan to display all 702 entries, then insert interface(rtablesize=[702,4])

 >
 >
 (1)
 >
 (2)
 >
 (3)
 >
 (4)
 >
 > Matrix([ lhs~(sol(0)),          seq(rhs~(sol(j)), j=0..7, 0.01)        ]       );
 (5)
 >

## You have multiple problems...

1. The dsolve() command accepts a "flat" list/set of equations and initial/boundary conditions - you have a set of three lists. This is a simple syntax error which can be corrected by writing dsolve({diffeqs[], initial_values[], final_conditions[]}, numeric, output = listprocedure);
2. Your system contains three first-order ODES, which means you can only have a grand total of three boundary/initial conditions. You have five (maybe 6, see next point)
3. You cannot have an inequality as a meaningful boundary condition,  so 0 < B(5.0) is completely unacceptable
4. Your system cannot include an entirely arbitrary function 'w(t)', with no derivatives or boundary/initial conditions

## Use the display() command...

as in the "toy" example attached.

 > restart;   with(plots):   with(plottools):   p1:=cylinder([0,0,0], 1, Pi, transparency=0.5, color=blue):   p2:=spacecurve([cos(2*t),sin(2*t),t],t=0..Pi, color=red, thickness=5):   display([p1,p2]);
 >

is that in the defintion of the function 'f', you have used both curly braces (ie {}) and square braces (ie []) for the simple grouping of terms. In Maple anything within {} designates a set and anythin within [] designates a member of an indexable quantity. For simple grouping of terms, use () only

## MAny ways to do this...

depending on the precise nature of the file to be imported. The attached worksheet will read a file containing

```[[0.1, 0.0769540597, 0.1477783335, 0.1393069312,

0.0763361154, 0.1477867626, 0.1393072151, 0.0763361266,

0.1477867830, 0.1393071934], [0.3, 0.1093424148, 0.1120401102,

0.1509302274, 0.1072278404, 0.1121142033, 0.1509369166,

0.1072278479, 0.1121142168, 0.1509369024], [0.5, 0.1392030568,

0.0853083077, 0.1558066181, 0.1353291378, 0.0855066806,

0.1558355785, 0.1353291558, 0.08550671332, 0.1558355439], [0.7,

0.1662374563, 0.0652194693, 0.1562235596, 0.1604342222,

0.0655908735, 0.1562974878, 0.1604342352, 0.06559089617,

0.1562974637], [0.9, 0.1903821623, 0.0500537619, 0.1537887672,

0.1825594352, 0.0506356391, 0.1539346707, 0.1825594528,

0.05063567192, 0.1539346352], [1.1, 0.2117168860, 0.0385555127,

0.1496220815, 0.2018541863, 0.0393727175, 0.1498707561,

0.2018542024, 0.03937274753, 0.1498707234], [1.3, 0.2303874000,

0.0298096396, 0.1444864012, 0.2185409755, 0.0308687065,

0.1448804021, 0.2185409880, 0.03086872986, 0.1448803765], [1.5,

0.2465077820, 0.0231661161, 0.1388614678, 0.2328759081,

0.0244336214, 0.1394893808, 0.2328759200, 0.02443364533,

0.1394893543]]```

and produce a matrix. You will have to change the 'pathname' and 'fileName' assignments to something appropriate for your installation.

 > restart;   pathName:="C:/Users/TomLeslie/Desktop/":   fileName:="dat.txt":   Data:=readline(StringTools:-Join( [pathName, fileName],"")):   line:=1:   while line<>0 do         line:=readline(StringTools:-Join( [pathName, fileName],"")):         if   `and`(length(line)>0, line<>0)         then Data:=StringTools:-Join( [Data, line]);         fi;   od:   myData:=Matrix(parse~(Data));
 (1)
 >
 >

## Some problems you are going to have...

1. Your pdes use the parameter 'phi', which is nowhere assigned a value, although the list of parameters contains two values for the parameter 'delta'. Maybe one of these is supposed to be a phi-value?
2. The parameter list has two values for 'nu'. Which ine do you want to use?
3. The supplied graphs are two-dimensional. Any graph of u(y,t) versus 'y' and 't' is obviously three-dimensional. Hence the 2D graphs you supply represent a "slice" through the 3D graph at a specific value of 't',. This value is nowhere given

When you resolve the above problema, the code you need will be something like that shown in the attached, where I have guessed that delta=0.1, phi=0.5, nu=1 and plotted u(y, 2). (Obviously these guesses may well be completely wrong!). The attached only supplies a graph for the variation of the parameter 'Pr'. Variation with respect to other parametrs would be performed similarly.

 > restart;   inf:=10:   pdes:= diff(u(y,t),t)-xi*diff(u(y,t),y)=diff(u(y,t),y\$2)/(1+lambda__t)+Gr*theta(y,t)+Gc*C(y,t)-M*u(y,t)-K*u(y,t),          diff(theta(y,t),t)-xi*diff(theta(y,t),y)=1/Pr*diff(theta(y,t),y\$2)+phi*theta(y,t),          diff(C(y,t),t)-xi*diff(C(y,t),y)=1/Sc*diff(C(y,t),y\$2)-delta*C(y,t)+nu*theta(y,t):   conds:= u(y,0)=0, theta(y,0)=0, C(y,0)=0,           u(0,t)=0, D[1](theta)(0,t)=-1, D[1](C)(0,t)=-1,           u(inf,t)=0, theta(inf,t)=0, C(inf,t)=0:   pars:= { Gr=1, Gc=1, M=1, nu=1, lambda__t=0.5,            Sc=0.78, delta=0.1, phi=0.5, K=0.5, xi=0.5          }
 (1)
 > PrVals:=[0.71, 1.00, 3.00, 7.00]:   colors:=[red, green, blue, black]:   for j from 1 to numelems(PrVals) do       pars1:=`union`( pars, {Pr=PrVals[j]}):       pdSol:= pdsolve( eval([pdes], pars1),                        eval([conds], pars1),                        numeric                      );       plt[j]:=pdSol:-plot( u(y,t), t=2, y=0..inf, numpoints=200, color=colors[j]);   od:   plots:-display( [seq(plt[j], j=1..numelems(PrVals))]);
 >
 >

## Well...

1. you use the name gamma in your ODEs, then define a parameter called localgamma - if these are supposed to be the same, then I suggesst using the latter name throughout
2. you use the name ic4 in one of your initial conditions, then define a parameter called ics4 - if these are supposed to be the same, then be consistent
3. the parameter N is used in your ODEs and is nowhere defined - it will have to be

In the attached, I have guessed your intent and fixed all of the above. It is up to you to determine whether my guesses are correct before further progress can be made.

 > restart; with(plots): with(plottools): with(DEtools): eqn1 := diff(V(t), t) = pi*p - (alpha + mu)*V(t),         V(0) = ic1; eqn2 := diff(S(t), t) = alpha*V(t) + (1 - p)*pi - beta*S(t)*In(t)/N - mu*S(t),         S(0) = ic2; eqn3 := diff(In(t), t) = beta*S(t)*In(t)/N - (mu + delta + localgamma)*In(t),         In(0) = ic3; eqn4 := diff(R(t), t) = localgamma*In(t) - mu*R(t),         R(0) = ic4; pi := 487845; p := 0.948; alpha := 0.054; beta := 0.955; mu := 0.005; delta := 0.03; localgamma := 0.935; ic1 := 484465; ic2 := 31999760; ic3 := 26305; ic4 := 12470; N:=10000; dsol := dsolve([eqn1, eqn2, eqn3, eqn4], numeric); odeplot( dsol,          [ [t, V(t), color = plum],            [t, S(t), color = blue],            [t, In(t), color = cyan],            [t, R(t), color = green]          ],          t = 0 .. 200,          thickness = 3        );
 >

## So far as I can tell...

Maple is producing two solutions, each of which is a two-element arctan() function, ie arctan(x,y).

Suggest arctan() help page for further enlightment

## SImilar problem to your last question...

Your wotksheet executes without error in Maple2023, but produces errors in your version (Maple2018). I tried (briefly) to make it function in Maple2018, but I failed.
Best suggestion I can come up with is not to use units at all, then use the labels option of the plot() command, to produce axes labels - as in the attached.

 >
 (1)
 >
 (2)
 >
 >
 >

## You can use maximize()...

provided that your expression is dimensionally correct. Nowhere in your original worksheeet did you specify that the variable 't' has 'seconds' as a unit. Even with this fix, the lost connection to kernel occurs in all versions of Maple up to Maple2022

See the  attached (in Maple2023 (which probably isn't very useful for you!)

 >
 >
 (1)
 >

## Observation...

Everything "works" (ie all plots are produced)  in Maple 2023, except for

```plot(f, t=0.9..1.12, adaptive=true);  #no graph
```

Either of the following "works" in Maple 2023

```plot(f, t=0.9..1.12, adaptive=geometric);  #no graph
```

## Not sure why...

yould wrie your own numeric method, for a differential equation which has an analytic solution (unless you want it for comparison purposes). The attached has twow figures

2. plot of the analytic solution
 >
 >
 >
 >
 (1)
 >
 (2)
 >
 (3)
 >
 (4)
 >
 >
 >

## In other words...

you want the coeffiicients of a polynomial? As in the attached.

 > restart;   f:= n->expand( (3+2*x+3*x^2)^n ); # # The polynomial #   f(7); # # Its coefficients #   coeffs( f(7), x);
 (1)
 >

## Were you hoping for something...

a bit more like the attached? (Arrays of plots do not render on this website - hence the download link only)

`eval(SensibilidadeK*(Profundidade, t, k, C), [k = 1, C = 1])`