Advanced Search

Maple Questions and Posts

These are Posts and Questions associated with the product, Maple
View unanswered posts
Maple Questions and Posts Feed

Accuracy issues with s-stage Runge-Kutta method...

Asked by: John_824 5
ode homework numeric
October 09 2023
1  2

Hi!

I have M number of linear differential equations. I have solved this system using the 1,2,3,4 stag RK method in the attached file but did not find a significant difference in the accuracy. Kindly see what's wrong there. 

Thank you!

s-stage.mw

preview print function not working correct after ...

Asked by: 2cUniverse 115
print export osx
October 09 2023
2  0

When printing a Maple Worksheet  often I go in the PrintPrewiev of the Mac and then select some sides to print.

This is not working anymore since I have updated from Maple 2021 to Maple 2023.

Is this known ?

Any help for his ?

Invalid left hand side in assignment - issue with ...

Asked by: NeraSnow 55
syntax assignment
October 08 2023
2  7

I was trying to assign a Vector to a Vector inside a procedure. However, when the Vector has a size of 5, I was able to do that. But not when the Vector has a size of 7.

I have the following example. I can't assign the Vector over when the Vector has a size of 7 as shown in oneStep_egcd2.  But I can do that when the size is 5 in oneStep_egcd. Copying the contents over using a for loop works as expected. 

> 

 

> 

a:= 17; b:= 5;

17

  

5

 (1)

 

> 

# <s, a, t, b, g>
prev := <1, a, 0, b, a>;

Vector(5, {(1) = 1, (2) = 17, (3) = 0, (4) = 5, (5) = 17})

 (2)
> 

curr := <0, a, 1, b, b>;

Vector(5, {(1) = 0, (2) = 17, (3) = 1, (4) = 5, (5) = 5})

 (3)
> 

oneStep_egcd := proc(prev::Vector, curr::Vector)
    local q, t1, t2, t3;
    if (curr[5] <> 0) then
        q := iquo(prev[5], curr[5]);
        t1, t2, t3 := prev[1] - curr[1] * q, prev[3] - curr[3] * q, prev[5] - curr[5] * q;
        prev = curr;
        curr[1], curr[3], curr[5] := t1, t2, t3;
    end if:
end proc:
> 

oneStep_egcd(prev, curr)

1, -3, 2

 (4)

NULL

NULL

> 

# <s, expr_s, a, t, expr_t, b, g>
prev := <1, 1, a, 0, 0, b, a>;

Vector(7, {(1) = 1, (2) = 1, (3) = 17, (4) = 0, (5) = 0, (6) = 5, (7) = 17})

 (5)
> 

curr := <0, 0, a, 1, 1, b, b>;

Vector(7, {(1) = 0, (2) = 0, (3) = 17, (4) = 1, (5) = 1, (6) = 5, (7) = 5})

 (6)
> 

oneStep_egcd2 := proc(prev::Vector, curr::Vector)
    local q, sb_q, t1, t2, t4, t5, t7;
    if (curr[7] <> 0) then
        q := iquo(prev[7], curr[7]);
        sb_q := (q);
        t1 := prev[1] - curr[1] * q;
        t2 := prev[2] - curr[2] * q;
        t4 := prev[4] - curr[4] * q;
        t5 := prev[5] - curr[5] * q;
        t7 := prev[7] - curr[7] * q;
        prev := curr;
        curr[1], curr[2], curr[4], curr[5], curr[7] := t1, t2, t4, t5, t7;
    end if:
end proc:
> 

oneStep_egcd2(prev, curr)

Error, (in oneStep_egcd2) invalid left hand side in assignment

  

NULL

> 

oneStep_egcd3 := proc(prev::Vector, curr::Vector)
    local q, sb_q, t1, t2, t4, t5, t7, i;
    if (curr[7] <> 0) then
        q := iquo(prev[7], curr[7]);
        sb_q := (q);
        t1 := prev[1] - curr[1] * q;
        t2 := prev[2] - curr[2] * q;
        t4 := prev[4] - curr[4] * q;
        t5 := prev[5] - curr[5] * q;
        t7 := prev[7] - curr[7] * q;
        for i to 7 do
            prev[i] := curr[i];
        end do;
        curr[1], curr[2], curr[4], curr[5], curr[7] := t1, t2, t4, t5, t7;
    end if:
end proc:
> 

oneStep_egcd3(prev, curr)

1, 1, -3, -3, 2

 (7)

NULL

Download ErrorExample.mw

Polynomial coefficients program code???...

Asked by: delvin 60
coefficients
October 08 2023
1  5

Hello

I want to write a program to get unknown coefficients of multiple polynomials. I have a problem with this program. The code sometimes doesn't work. Can anyone help me? It's very important to me.

> 

restart
> 

with(student)

``

EQ[0] := tanh(d)*b[1]*(b[1]+1)

tanh(d)*b[1]*(b[1]+1)

 (1)

EQ[1] := -(-1+(a[1]-b[1]-1)*tanh(d)^2+(a[0]+1)*tanh(d))*b[1]

-(-1+(a[1]-b[1]-1)*tanh(d)^2+(a[0]+1)*tanh(d))*b[1]

 (2)

EQ[2] := tanh(d)*((a[1]-b[1])*(a[0]+1)*tanh(d)-b[1]^2-a[1])

tanh(d)*((a[1]-b[1])*(a[0]+1)*tanh(d)-b[1]^2-a[1])

 (3)

EQ[3] := (-a[1]^2+(2*b[1]-1)*a[1]-b[1]^2-b[1])*tanh(d)^2+(a[1]+b[1])*(a[0]+1)*tanh(d)-a[1]-b[1]

(-a[1]^2+(2*b[1]-1)*a[1]-b[1]^2-b[1])*tanh(d)^2+(a[1]+b[1])*(a[0]+1)*tanh(d)-a[1]-b[1]

 (4)

EQ[4] := -tanh(d)*((a[1]-b[1])*(a[0]+1)*tanh(d)+a[1]^2+b[1])

-tanh(d)*((a[1]-b[1])*(a[0]+1)*tanh(d)+a[1]^2+b[1])

 (5)

EQ[5] := -a[1]*(-1+(-a[1]+b[1]-1)*tanh(d)^2+(a[0]+1)*tanh(d))

-a[1]*(-1+(-a[1]+b[1]-1)*tanh(d)^2+(a[0]+1)*tanh(d))

 (6)

EQ[6] := (a[1]+1)*a[1]*tanh(d)

(a[1]+1)*a[1]*tanh(d)

 (7)

Eqs := {seq(EQ[i], i = 0 .. 6)}

Sol := solve(Eqs, {a[0], a[1], b[1]})NULL

{a[0] = a[0], a[1] = 0, b[1] = 0}, {a[0] = (tanh(d)^2-tanh(d)+1)/tanh(d), a[1] = -1, b[1] = -1}, {a[0] = -(tanh(d)-1)/tanh(d), a[1] = -1, b[1] = 0}, {a[0] = -(tanh(d)-1)/tanh(d), a[1] = 0, b[1] = -1}

 (8)

for i from 2 to 4 do Case[i] := allvalues(Sol[i]) end do

{a[0] = -(tanh(d)-1)/tanh(d), a[1] = 0, b[1] = -1}

 (9)

``

T1.mw

T2.mw

Reading a PNG (or SVG)-file, adding a white frame ...

Asked by: 2cUniverse 115
image export imagetools
October 08 2023
2  11

My Maple-Program is creating svg-files (or PNG) in a squared format.
There are a lot of these files in a local folder.
I would like to add a white surrounding frame (10% of the Imagesize) to it and afterwards store it as a TIFF file.

Studying ImageTools I have not found a way to solve this with ImageTools.

Until now I am using a graphic software to add a frame and store as TIFF in a desired resolution.

Now I plan to do this process automatically with maple.

Thanks for support :)

Integro-Diff-Eq steps to solution...

Asked by: Joseph Poveromo 25
integro-differential-equation
October 08 2023
0  2

For the following Equation:

Equation := int(diff(u(x), x)*v(x), x) = int(u(x)^(1/2)*v(x), x)^(-2/3);
Maplesoft finds the following solution:

Solution1:=3/4*u(x)^(4/3) + 2/3*u(x)^(5/6)*Intat(1/(sqrt(u(x))*Int(v(_b), _b))^(5/3), _b = x) + _C1 = 0

or , which I believe as an alternative, can be written as

Solution2:=3/4*u(x)^(4/3) + 2/3*u(x)^(5/6)*Int(1/(sqrt(u(x))*Int(v(x),x))^(5/3) +_C1=0

My question is how did Maple arrive at 'Solution1' from 'Equation'? In other words, can someone fill

in the steps between 'Equation'  and 'Solution1'? Or even, prove that Solution 1 is a valid solution to Equation.

Plugging the Solution1 into Equation, did not clearly demonstate the validity of the solution (to me at least)

Unfortunately, I am still unable to post the corresponding Maplesoft worksheet onto this post.

Gfun's listtoalgeq erroneous output?...

Asked by: AlexShura 25
gfun oeis
October 07 2023
0  5

Invoked by the OEIS superseeker, Maple "gfun" package "listtoalgeq" identified possible lgdegf for https://oeis.org/A035001

1, 2, 4, 5, 8, 12, 14, 16, 28, 32, 37, 64, 94, 106, 128, 144, 232, 256, 289, 320, 512, 560, 704, 760, 838, 1024, 1328, 1536, 1944, 2048, 2329, 3104, 3328, 4096, 4864, 6266, 6802, 7168, 8192, 11952, 15360, 16384, 16428, 19149, 28928, 32768, 37120, 42168 

as follows:

1024-5120*a(n)+11520*a(n)^2-15360*a(n)^3+13440*a(n)^4-8064*a(n)^5+3360*a(n)^6-960*a(n)^7+180*a(n)^8-20*a(n)^9+a(n)^{10}

The coefficients of above polynomial are:

{1, -20, 180, -960, 3360, -8064,13440, -15360, 11520, -5120, 1024,...}

It is interesting that the absolute values of above polynomial coefficients satisfy a(n) of

https://oeis.org/A013609

for n=55...65,

which is the 11th row in the triangle presentation of A013609, so in other words the absolute values of above polynomial coefficients are T={11, k} for k=1...11

System of equation into matrix for each "n"...

Asked by: Muhammad Usman 240
matrix
October 07 2023
0  2

Dear Users!

I hope you are doing well. I have the following discretized form

for n>=1 and j=0..M. We obtained the following matrix equation for any "n" and j=0..M as:

I want matrix proc of any useful way to define A^n, u^n, and b^n. I am waiting for your positive response. Thanks in advancs

How to create geom3d point functions?...

Asked by: rowlesmr 50
geom3d point
October 07 2023
1  1

I would like to do some point and line calculations, and the geom3d package seems a good one.

But how do I define a generic point such that I can use it in a function to generate genric points and lines?

eg

    P:=(a,b,c)->point(P, a,b,c) 

or equivalent?

I would then want to define generic lines between generic points.

Maple gives trivial solution to linear system with...

Asked by: zenterix 405
linear-algebra parametric
October 06 2023
1  1

Consider the following. 

We have a matrix A with five parameters in it.

If we tell Maple to solve Ax=0, it gives us the trivial solution.

Then, if I make a matrix B where I choose specific values for those five parameters and I tell Maple to solve Bx=0 I get a non-trivial solution.

Why doesn't Maple give me a more informative result in the Ax=0 case? 

A := Matrix(5, 5, {(1, 1) = 1, (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (1, 5) = a__1, (2, 1) = 0, (2, 2) = 1, (2, 3) = 0, (2, 4) = 0, (2, 5) = b__1, (3, 1) = -216, (3, 2) = -18, (3, 3) = 0, (3, 4) = 0, (3, 5) = c__1, (4, 1) = 0, (4, 2) = 0, (4, 3) = 1, (4, 4) = 0, (4, 5) = d__1, (5, 1) = 0, (5, 2) = 0, (5, 3) = 0, (5, 4) = 1, (5, 5) = e__1}) = Matrix(%id = 36893488151959194068)NULL

with(LinearAlgebra)

GaussianElimination(A)

Matrix(%id = 36893488151959184316)

 (1)

LinearSolve(A, `<,>`(0, 0, 0, 0, 0))

Vector[column](%id = 36893488151959183100)

 (2)

B := Matrix(5, 5, {(1, 1) = 1, (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (1, 5) = 0, (2, 1) = 0, (2, 2) = 1, (2, 3) = 0, (2, 4) = 0, (2, 5) = 0, (3, 1) = -216, (3, 2) = -18, (3, 3) = 0, (3, 4) = 0, (3, 5) = 0, (4, 1) = 0, (4, 2) = 0, (4, 3) = 1, (4, 4) = 0, (4, 5) = 1, (5, 1) = 0, (5, 2) = 0, (5, 3) = 0, (5, 4) = 1, (5, 5) = 4}) = Matrix([[1, 0, 0, 0, 0], [0, 1, 0, 0, 0], [-216, -18, 0, 0, 0], [0, 0, 1, 0, 1], [0, 0, 0, 1, 4]])NULL

GaussianElimination(B) = Matrix(%id = 36893488151959174788)NULL

LinearSolve(B, `<,>`(0, 0, 0, 0, 0))

Vector[column](%id = 36893488151959164916)

 (3)

``

NULL

An even easier way to show a non-trivial solution is to just make the entire last column zero.

Download polyn_matrix.mw

simplify piecewise function ...

Asked by: Zeineb 540
simplify piecewise incomplete-question
October 06 2023
0  2

Dear all 

I have a function defined on many sub-intervals, how can I simplify this the funciton obtained at each iteration. I hope obatin B_{i,1}, B_{i,2}, and B_{i,3} 

B_Spline.mw

Thank you 

Extract vectors from a vector...

Asked by: Muhammad Usman 240
vector
October 06 2023
0  2

Dear Users!

I hope everyone is fine here. I have three vectors say V[1], V[2] and V[3] as:

restart; with(LinearAlgebra); with(linalg);
V[1] := Vector([3, 2, -4]);
V[2] := Vector([1, 2/3, 7]);
V[3] := Vector([-9, 0, 1/2]);

I can define a Vector V have V[1], V[2] and V[3] using blockmatrix command as:

C := blockmatrix(3, 1, [V[1], V[2], V[3]]);

My question is that how I can extract vectors V[1], V[2] and V[3] from C? Thanks in advance

Four-vectors in Maple...

Asked by: Traruh Synred 55
vector physics duplicate-question
October 05 2023
1  19

Is there a simple 4-vector package available for Maple? Yes, I know about using Tensors, but don't want to fool with metrics and raising and lowering operators.

I have written my  own, but it seems buggy. I generate some vectors in in CM that should add up to a massless neutrino, but some times when I boost them to the 'lab' and do the addition there the neutrion acquires a substantial mass. It's not a merely numerical issue. It's weird because most of the time the neutrion mass is 0 or close to it.

My homemade boost proc 'boostTo' is attached.

Error, (in minimize) unexpected options: {x}...

Asked by: Missis 15
extrema minimize
October 05 2023
1  5

Привет. У меня проблема с пестремой, помогите пожалуйста:

restart;
y:=sqrt((x^2*(x^2-4))/(x^2-1));
readlib(extrema): readlib(maximize): readlib(minimize):
extrema(y,{},x,'s');s;
ymax:=maximize(y,{x});
Error, (in maximization) unexpected parameters: {x}
ymin:=minimize(y,{x});
Error, (in minimization) unexpected parameters: {x}

unable to plot 3d and 2d plots...

Asked by: KIRAN SAJJAN 40
ode numeric differential-equations
October 05 2023
0  3

  I am unable to draw both 3d plots sowing error please help me to solve

> 

restart:NULLNULL
> 

p1 := 0.1e-1; p2 := 0.2e-1; p3 := 0.1e-1; Px := p1+p2+p3
> 

rf := 1050; kf := .52; cpf := 3617; sigmaf := .8
> 

sigma1 := 25000; rs1 := 5200; ks1 := 6; cps1 := 670
> 

sigma2 := 59.7*10^6; rs2 := 8933; ks2 := 400; cps2 := 385
> 

sigma3 := 2380000; rs3 := 4250; ks3 := 8.9538; cps3 := 686.2
> 

NULL
> 

B1 := 1+2.5*Px+6.2*Px^2; B2 := 1+13.5*Px+904.4*Px^2; B3 := 1+37.1*Px+612.6*Px^2; B4 := (ks1+2*kf-2*Px*(kf-ks1))/(ks1+2*kf+Px*(kf-ks1)); B5 := (ks2+3.9*kf-3.9*Px*(kf-ks2))/(ks2+3.9*kf+Px*(kf-ks2)); B6 := (ks3+4.7*kf-4.7*Px*(kf-ks3))/(ks3+4.7*kf+Px*(kf-ks3))
> 

a2 := B1*p1+B2*p2+B3*p3
> 

a1 := 1-p1-p2-p3+p1*rs1/rf+p2*rs2/rf+p3*rs3/rf
> 

a3 := 1-p1-p2-p3+p1*rs1*cps1/(rf*cpf)+p2*rs2*cps2/(rf*cpf)+p3*rs3*cps3/(rf*cpf)
> 

a4 := B4*p1+B5*p2+B6*p3
> 

NULL
> 

a5 := 1+3*((p1*sigma1+p2*sigma2+p3*sigma3)/sigmaf-p1-p2-p3)/(2+(p1*sigma1+p2*sigma2+p3*sigma3)/((p1+p2+p3)*sigmaf)-((p1*sigma1+p2*sigma2+p3*sigma3)/sigmaf-p1-p2-p3))

``

``

> 



NULL
> 

ODE:=[(a2+K)*(diff(U0(eta), eta, eta))/a1-Ra*(diff(U0(eta), eta))+lambda0/a1-a5*M1^2*U0(eta)/a1+K*(diff(N0(eta), eta))/a1+la*Ra*Theta0(eta)*(1+Qc*Theta0(eta)), (a2+K)*(diff(U1(eta), eta, eta))/a1-H^2*l1*U1(eta)-Ra*(diff(U1(eta), eta))+lambda1/a1-a5*M1^2*U1(eta)/a1+K*(diff(N1(eta), eta))/a1+la*Ra*(Theta1(eta))(1+2*Qc*Theta0(eta)), diff(N0(eta), eta, eta)-Ra*a1*Pj*(diff(N0(eta), eta))-2*n1*N0(eta)-n1*(diff(U0(eta), eta)), diff(N1(eta), eta, eta)-Ra*a1*Pj*(diff(N1(eta), eta))-2*n1*N1(eta)-n1*(diff(U1(eta), eta))-H^2*a1*Pj*l1*N1(eta), (a4/(a3*Pr)-delta*Ra^2/H^2+4*Rd*(1+(Tp-1)^3*Theta0(eta)^3+3*(Tp-1)^2*Theta0(eta)^2+(3*(Tp-1))*Theta0(eta))/(3*a3*Pr))*(diff(Theta0(eta), eta, eta))-Ra*(diff(Theta0(eta), eta))+a5*Ec*M1^2*U0(eta)^2/a3+(a2+K)*Ec*(diff(U0(eta), eta))^2/a1+Q*Theta0(eta)/a3+4*(diff(Theta0(eta), eta))^2*Rd*(3*(Tp-1)+6*(Tp-1)^2*Theta0(eta)+3*(Tp-1)^3*Theta0(eta)^2)/(3*a3*Pr), (a4/(a3*Pr)-delta*Ra^2/H^2+4*Rd*(1+(Tp-1)^3*Theta0(eta)^3+3*(Tp-1)^2*Theta0(eta)^2+(3*(Tp-1))*Theta0(eta))/(3*a3*Pr))*(diff(Theta1(eta), eta, eta))-(H^2*l1+2*Ra*delta*l1+Ra)*(diff(Theta1(eta), eta))+(Q/a3-delta*H^2*l1^2)*Theta1(eta)+2*(a2+K)*Ec*(diff(U0(eta), eta))*(diff(U1(eta), eta))/a1+2*a5*Ec*M^2*U0(eta)*U1(eta)/a3+4*(diff(Theta0(eta), eta, eta))*Theta1(eta)*Rd*(3*(Tp-1)+6*(Tp-1)^2*Theta0(eta)+3*(Tp-1)^3*Theta0(eta)^2)/(3*a3*Pr)+4*Rd*(diff(Theta0(eta), eta))^2*(6*(Tp-1)^2*Theta1(eta)+6*(Tp-1)^3*Theta0(eta)*Theta1(eta))/(3*a3*Pr)+4*Rd*(diff(Theta1(eta), eta))*(diff(Theta0(eta), eta))*(6*(Tp-1)+6*(Tp-1)^3*Theta0(eta)^2+12*(Tp-1)^2*Theta0(eta))/(3*a3*Pr)]:
> 


(LB,UB):= (0,1):


BCs:= [
  
  U0(0) = 0, U1(0) = 0, N0(0) = 0, N1(0) = 0, Theta0(0) = 0, Theta1(0) = 0, U0(1) = 0, U1(1) = 0, N0(1) = 0, N1(1) = 0, Theta0(1) = 1, Theta1(1) = 0
]:
> 

NULL
> 


Params:= Record(
   
   M1=  1.2, Rd=0.8,la=0.8,n1=1.2,Q=0.2,Pj=0.001,Ra=0.8,Ec=1,    Pr= 21,   delta= 0.2,    t1= (1/4)*Pi, lambda0=2,lambda1=3,   Qc= 0.1,    l1= 1,K=0.4,H=3 ,deltat=0.05  ):
   
> 

NBVs:= [   
 
a1**D(U0)(0) = `C*__f` , # Skin friction coefficient
 (a4+(4*Rd*(1/3))*(1+(Tp-1)*(Theta0(0)+0.1e-2*exp(l1*t1)*Theta1(0)))^3)*((D(Theta0))(0)+0.1e-2*exp(l1*t1)*(D(Theta1))(0)) = `Nu*`    # Nusselt number     
]:
Nu:= `Nu*`:
Cf:= `C*__f`:

 
> 

Solve:= module()
local
   nbvs_rhs:= rhs~(:-NBVs), #just the names
   Sol, #numeric dsolve BVP solution of any 'output' form
   ModuleApply:= subs(
      _Sys= {:-ODEs[], :-BCs[], :-NBVs[]},
      proc({
          M1::realcons:=  Params:-M1,
         Pr::realcons:= Params:-Pr,
         Rd::realcons:= Params:-Rd,
         la::realcons:= Params:-la,
         Tp::realcons:= Params:-Tp,
         n1::realcons:= Params:-n1,
         Q::realcons:= Params:-Q,
         Pj::realcons:= Params:-Pj,
         Ra::realcons:= Params:-Ra,
         Ec::realcons:= Params:-Ec,
         t1::realcons:=  Params:-t1,
         delta::realcons:= Params:-delta,
         lambda0::realcons:= Params:-lambda0,
         lambda1::realcons:= Params:-lambda1,
         Qc::realcons:= Params:-Qc,
         K::realcons:= Params:-K,
         l1::realcons:= Params:-l1,
         H::realcons:= Params:-H
      })
         Sol:= dsolve(_Sys, _rest, numeric);
         AccumData(Sol, {_options});
         Sol
      end proc
   ),
   AccumData:= proc(
      Sol::{Matrix, procedure, list({name, function}= procedure)},
      params::set(name= realcons)
   )
   local n, nbvs;
      if Sol::Matrix then
         nbvs:= seq(n = Sol[2,1][1,Pos(n)], n= nbvs_rhs)
      else
         nbvs:= (nbvs_rhs =~ eval(nbvs_rhs, Sol(:-LB)))[]
      fi;
      SavedData[params]:= Record[packed](params[], nbvs)
   end proc,
   ModuleLoad:= eval(Init);
export
   SavedData, #table of Records
   Pos, #Matrix column indices of nbvs
   Init:= proc()
      Pos:= proc(n::name) option remember; local p; member(n, Sol[1,1], 'p'); p end proc;
      SavedData:= table();
      return
   end proc ;
   ModuleLoad()
end module:
 
> 


 
> 

 

> 

#procedure that generates 3-D plots (dropped-shadow contour + surface) of an expression


ParamPlot3d:= proc(
   Z::{procedure, `module`}, #procedure that extracts z-value from Solve's dsolve solution
   X::name= range(realcons), #x-axis-parameter range
   Y::name= range(realcons), #y-axis-parameter range
   FP::list(name= realcons), #fixed values of other parameters
   {
      #fraction of empty space above and below plot (larger "below"
      #value improves view of dropped-shadow contourplot):
      zmargin::[realcons,realcons]:= [.05,0.15],
      eta::realcons:= :-LB, #independent variable value
      dsolveopts::list({name, name= anything}):= [],
      contouropts::list({name, name= anything}):= [],
      surfaceopts::list({name, name= anything}):=[]    
   }
)
local
   LX:= lhs(X), RX:= rhs(X), LY:= lhs(Y), RY:= rhs(Y),
   Zremember:= proc(x,y)
   option remember; #Used because 'grid' should be the same for both plots.
      Z(
         Solve(
            LX= x, LY= y, FP[],
            #Default dsolve options can be changed by setting 'dsolveopts':
            'abserr'= 0.5e-7, 'interpolant'= false, 'output'= Array([eta]),  
            dsolveopts[]
         )
      )
   end proc,
   plotspec:= (Zremember, RX, RY),
   C:= plots:-contourplot(
      plotspec,
      #These default plot options can be changed by setting 'contouropts':
      'grid'= [25,25], 'contours'= 5, 'filled',
      'coloring'= ['yellow', 'orange'], 'color'= 'green',
      contouropts[]
   ),
   P:= plot3d(
      plotspec,
      #These default plot options can be changed by setting 'surfaceopts':
      'grid'= [25,25], 'style'= 'surfacecontour', 'contours'= 6,
      surfaceopts[]
   ),
   U, L #z-axis endpoints after margin adjustment
;
   #Stretch z-axis to include margins:
   (U,L):= ((Um,Lm,M,m)-> (M*(Lm-1)+m*Um, M*Lm+m*(Um-1)) /~ (Um+Lm-1))(
      zmargin[],
      (max,min)(op(3, indets(P, 'specfunc'('GRID'))[])) #actual z-axis range
   );
   plots:-display(
      [
         plots:-spacecurve(
            {
               [[lhs(RX),rhs(RY),U],[rhs(RX),rhs(RY),U],[rhs(RX),rhs(RY),L]], #yz backwall
               [[rhs(RX),rhs(RY),U],[rhs(RX),lhs(RY),U],[rhs(RX),lhs(RY),L]]  #xz backwall
            },
            'color'= 'grey', 'thickness'= 0
         ),
         plottools:-transform((x,y)-> [x,y,L])(C), #dropped-shadow contours
         P
      ],
      #These default plot options can be changed simply by putting the option in the
      #ParamPlot3d call:
      'view'= ['DEFAULT', 'DEFAULT', L..U], 'orientation'= [-135, 75], 'axes'= 'frame',
      'labels'= [lhs(X), lhs(Y), Z], 'labelfont'= ['TIMES', 'BOLDOBLIQUE', 16],
      'caption'= nprintf(cat("%a = %4.2f, "$nops(FP)-1, "%a = %4.2f"), (lhs,rhs)~(FP)[]),
      'captionfont'= ['TIMES', 14],
      'projection'= 2/3,   
      _rest
   )
end proc:

> 

NULL
> 

NULL
> 

GetNu := proc (Sol::Matrix) options operator, arrow; Sol[2, 1][1, Solve:-Pos(:-Nu)] end proc
> 

ParamPlot3d(
   GetNu,Q= 0..5, Rd= 0..5, [
   
   Pr= 21   ],
   labels= [Q, gamma, Nu]
);

Error, (in plot/iplot2d/levelcurve) could not evaluate expression

  

``

Download P6_3D_plots.mw

First 163 164 165 166 167 168 169 Last Page 165 of 2217