Advanced Search

MaplePrimes Questions

How to use goto without “System error”...

Asked by: lcz 994
programming goto
June 04 2023
0  12

Let's put aside the drawbacks of using goto for now. Below is the program I wrote. When I run it, it seems to execute, but the Maple interface displays the message "System error,,at top level". Does this trigger an internal warning in Maple?

for i to 4 do
    for j to 4 do
         print([i,j]):
        if i = 2 and j = 3 then
             goto(al):
        end if:
    end do:
end do:
al;

                             [1, 1]

                             [1, 2]

                             [1, 3]

                             [1, 4]

                             [2, 1]

                             [2, 2]

                             [2, 3]

                  System error, , at top level

                               al

Using coeff and subs function...

Asked by: ijuptilk 100
substitution subs
June 03 2023
1  5

I have the following code below:

eq := -((2*(alpha[3]+alpha[2]))*((1/2)*cos(Theta(z, t))^2-(1/2)*sin(Theta(z, t))^2)+(alpha[3]-alpha[2])*(sin(Theta(z, t))^2+cos(Theta(z, t))^2))*(diff(U(z, t), z))-(alpha[3]-alpha[2])*(2*sin(Theta(z, t))^2*(diff(Theta(z, t), t))+2*cos(Theta(z, t))^2*(diff(Theta(z, t), t)))+(diff(Theta(z, t), z, z))*k

simplify(subs(simplify(coeff(eq, diff(U(z, t), z))) = t1, eq))

Please, why is maple not assigning t1 to the coefficient of diff(U(z, t), z)) in the expression? 

How to integrate this ode and how to substitue thi...

Asked by: C_R 1960
integration ode differential-equations
June 03 2023
1  9

The example below is an attempt to solve a boundary value problem (BVP) without manually replacing terms in expressions and some questions where there might be room for improvement. (Only the first integration step is shown to keep it short)

My objective is to have clean Maple input with minimal manual interventions or low-level manipulations (such as using op and similar commands).

NULL

> 

ode := diff(phi(s), s, s)+K*cos(phi(s)) = 0

diff(diff(phi(s), s), s)+K*cos(phi(s)) = 0

 (1)

NULL

> 

map(int, ode, s)

int(diff(diff(phi(s), s), s)+K*cos(phi(s)), s) = 0

 (2)
> 

Student[ODEs]:-Integrate(ode, phi(s))

Error, (in assuming) when calling 'Student:-ODEs:-Integrate'. Received: 'The given ODE, diff(diff(phi(s),s),s)+K*cos(phi(s)) = 0, is not exact'

  

...not a clue for me


Since I know the result: Multiplying by diff(phi(x), x) to make integration work.  Q1: Are there other ways to integrate a single step without a priori knowledge?

> 

(diff(phi(s), s))*ode

(diff(phi(s), s))*(diff(diff(phi(s), s), s)+K*cos(phi(s))) = 0

 (3)

NULL

> 

int(lhs((diff(phi(s), s))*(diff(diff(phi(s), s), s)+K*cos(phi(s))) = 0), s)+c__1 = 0

(1/2)*(diff(phi(s), s))^2+K*sin(phi(s))+c__1 = 0

 (4)

Boundary condition for s = L 

> 

Eval(diff(phi(s), s), s = L) = `φ'`[0]

Eval(diff(phi(s), s), s = L) = `φ'`[0]

 (5)

(4) evaluated at the boundary

> 

expand(map(Eval, (1/2)*(diff(phi(s), s))^2+K*sin(phi(s))+c__1 = 0, s = L))

Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)

 (6)

Q2: How to simplify expression (6) with (5) using Maple commands to get this expression for the integration constant

``

> 

c__1 = -(1/2)*`φ'`[0]^2-K*sin(Eval(phi(s), s = L))

c__1 = -(1/2)*`φ'`[0]^2-K*sin(Eval(phi(s), s = L))

 (7)

Some attempts

> 

algsubs(Eval(diff(phi(s), s), s = L) = `φ'`[0], Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)); applyrule(Eval(diff(phi(s), s), s = L) = `φ'`[0], Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)); simplify(Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L), {Eval(diff(phi(s), s), s = L) = `φ'`[0]}); Physics:-Substitute(Eval(diff(phi(s), s), s = L) = `φ'`[0], Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L))

Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)

  

Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)

  

Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)

  

Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)

 (8)

``

> 

map(`^`, Eval(diff(phi(s), s), s = L) = `φ'`[0], 2)

(Eval(diff(phi(s), s), s = L))^2 = `φ'`[0]^2

 (9)

Doing it manually

isolate(subs({Eval(0, s = L) = 0, Eval((1/2)*(diff(phi(s), s))^2, {s = L}) = (1/2)*`φ'`[0]^2}, Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L)), c__1)

c__1 = -(1/2)*`φ'`[0]^2-(Eval(K*sin(phi(s)), {s = L}))

 (10)

or

isolate(algsubs(value(Eval(diff(phi(s), s), s = L) = `φ'`[0]), value(Eval((1/2)*(diff(phi(s), s))^2, {s = L})+Eval(K*sin(phi(s)), {s = L})+c__1 = Eval(0, s = L))), c__1)

c__1 = -(1/2)*`φ'`[0]^2-K*sin(phi(L))

 (11)

subs(phi(L) = Eval(phi(s), s = L), c__1 = -(1/2)*`φ'`[0]^2-(Eval(K*sin(phi(s)), {s = L})))

c__1 = -(1/2)*`φ'`[0]^2-(Eval(K*sin(phi(s)), {s = L}))

 (12)

Q3: Are there ways to prevent commands evaluating Eval(phi(s), s = L)to phi(L)which might suggest incorrectly a new independent variable?

Q4: In (6), why is s = L sometimes in brackets {} ?

``
Download BVP_with_questions.mw

how do I solve this sum?...

Asked by: Harry Xu 20
sum homework do-my-homework-for-me
June 03 2023
0  2

how do I solve this one?...

Asked by: Harry Xu 20
homework differential-equations do-my-homework-for-me
June 03 2023
0  2

effect of evala on expression when using integrate...

Asked by: nm 8552
integration int evala
June 02 2023
1  6

Could someone help me understand why Maple's int hangs when using int on an integrand if evala is not first applied to the integrand?

When I mean hangs, I mean really hangs. I have timelimit and maple server.exe runs at high CPU for hrs. When adding evala(), it completes instantly and returns unevaluated integral (which is just fine as I did not expect this to be evaluated) but it does not hang which is the main point.

Is this considered normal behaviour or a bug?  Should then one always add evala() on the integrand before using int?  As I get many hangs on int (even when using timelimit). I am not familiar with evala command as I have never had to use it before, but by chance I tried it and noticed this.

V 2023 on windows 10. Worksheet attached.


 

186932

> 

restart;

186932

> 

expr:=-1/3*2^(2/3)/((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)*((-1/2*a^2*p-1/2*p^2*a+1/2*p*((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)*a-3)*(-p^2*(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p))^(2/3)+2^(2/3)*p*(-(a+3/2*p)*p*((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a^2*p+5/2*p^2*a+3/2*p^3+3))/(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p)/(-1/2*2^(2/3)*(-p^2*(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p))^(2/3)+p*(p*(-p^2*(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p))^(1/3)+2^(1/3)))/p^2:
try
    anti := timelimit(30,int(evala(expr),p));
catch:
    print("timed out OK");
end try:
print("done");

"done"

> 

restart;
> 

expr:=-1/3*2^(2/3)/((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)*((-1/2*a^2*p-1/2*p^2*a+1/2*p*((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)*a-3)*(-p^2*(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p))^(2/3)+2^(2/3)*p*(-(a+3/2*p)*p*((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a^2*p+5/2*p^2*a+3/2*p^3+3))/(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p)/(-1/2*2^(2/3)*(-p^2*(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p))^(2/3)+p*(p*(-p^2*(-((a^2*p+2*a*p^2+p^3+4)/p)^(1/2)+a+p))^(1/3)+2^(1/3)))/p^2:
try
    anti := timelimit(30,int(expr,p));
catch:
    print("timed out");
end try;
print("done");


The second case above just hangs. Make sure to save all your work. I found I have to terminate all of Maple sometimes as killing server.exe seems to make the frontend hangs and not respond, ending up losing work in other worksheets which I did not save.

Download int_hangs_june_2_2023.mw

Error, (in dsolve/numeric/bvp/convertsys) unable t...

Asked by: Saha 70
ode homework dsolve
June 02 2023
0  2

I'm new in maple.I'm just learning, i couldn't get the solution.Please help me to solve this Error.

Here is my code.

sps-RK_Method.mw

Why is the output of an inert integral not in 2D M...

Asked by: C_R 1960
integral inert 2d-output
June 02 2023
1  3

I am updating older files.

Executing with the attached code snippet from a worksheet created with Maple 16 outputs for an inert integral

instead of

Why is that and how can I fix my old worksheets to make them work with Maple 2023?

Maybe related: Execution with (or return) does not evaluate the document blocks. When all document blocks are expanded with "right click show command" the cursor does not advance to the next execution group. I can't remenber if this behaviour was allways like this or has changed.

A suggestion for Maplesoft: I have stopped using document blocks with hidden code since hiding and showing commands requires too much time (to many clicks and mouse movements). A simple double click on the marker of the document block could facilitate hiding and showing code allot.

Output_of_Int_not_in_2d.mw

Re: Substituting an expression...

Asked by: pik1432 245
substitution
June 02 2023
2  8

Hello there, 

Is there any chance to ask for your expertise in this issue of substitution?

I tried to make the expression 'eq_m1' to 'eq_m1_desired' by substituting 'L__ads*L__fd/(L__fd + L__ads)' as 'L__ads_p', but did not succeeded, after a series of different attempts. Would you show me how to achive the desired expression?

> 

restart;
> 

with(LinearAlgebra):
> 

with(DynamicSystems):
> 

interface(imaginaryunit=j):
> 

eq_m1 := m1 = (X__Tq*E__B*sin(delta) - R__T*E__B*cos(delta))/(R__E^2 + 2*R__a*R__E + R__a^2 + X__E^2 + X__E*L__ads*L__fd/(L__fd + L__ads) + 2*X__E*L__l + L__aqs*X__E + L__aqs*L__ads*L__fd/(L__fd + L__ads) + L__aqs*L__l + L__l*L__ads*L__fd/(L__fd + L__ads) + L__l^2);

m1 = (X__Tq*E__B*sin(delta)-R__T*E__B*cos(delta))/(R__E^2+2*R__a*R__E+R__a^2+X__E^2+X__E*L__ads*L__fd/(L__fd+L__ads)+2*X__E*L__l+L__aqs*X__E+L__aqs*L__ads*L__fd/(L__fd+L__ads)+L__aqs*L__l+L__l*L__ads*L__fd/(L__fd+L__ads)+L__l^2)

 (1)
> 

eq_m1_desired := m1 = (X__Tq*E__B*sin(delta) - R__T*E__B*cos(delta))/(R__E^2 + 2*R__a*R__E + R__a^2 + X__E^2 + X__E*L__ads_p + 2*X__E*L__l + L__aqs*X__E + L__aqs*L__ads_p + L__aqs*L__l + L__l*L__ads_p + L__l^2);

m1 = (X__Tq*E__B*sin(delta)-R__T*E__B*cos(delta))/(L__aqs*L__l+L__aqs*X__E+L__aqs*L__ads_p+L__l^2+2*L__l*X__E+L__l*L__ads_p+R__E^2+2*R__E*R__a+R__a^2+X__E^2+X__E*L__ads_p)

 (2)
> 

######
> 

#eq_m1a := subs({L__ads*L__fd/(L__fd + L__ads) = L__ads_p}, eq_m1);
> 

#eq_m1a := subs({L__ads_p = L__ads*L__fd/(L__fd + L__ads)}, eq_m1);
> 

#eq_m1a := applyrule(L__ads*L__fd/(L__fd + L__ads) = L__ads_p, eq_m1);
> 

#eq_m1a := algsubs(L__ads_p = (L__ads*L__fd)/(L__fd + L__ads), eq_m1);
> 

#eq_m1a := algsubs((L__ads*L__fd)/(L__fd + L__ads) = L__ads_p, eq_m1);
> 

#eq_m1a := eval(eq_m1, L__ads_p = (L__ads*L__fd)/(L__fd + L__ads));
> 

 

Download Q20230602.mw

On the bondage number of a graph....

Asked by: lcz 994
graphtheory
June 02 2023
2  4
  • A set D of vertices in a graph G is a dominating set if each vertex of G that is not in D is adjacent to at least one vertex of D.
  • The minimum cardinality among all dominating sets in G is called the domination number of G and denoted σ(G).
  • We define the bondage number b(G) of a graph G to be the cardinality of a smallest set E of edges for which σ(G−E)>σ(G). 

In the previous post, Carl Love provided a code for calculating domination number, so we could use it to calculate bondage number. I wrote a piece of code using a relatively basic approach. I feel like there is room for improvement.  I still can't determine the bondage number of the following graph so far.

 
> 

restart:
with(GraphTheory):
with(Iterator):
#Author: Carl Love <carl.j.love@gmail.com>, the Glorious 1st of June, 2023
#
Domination_number:= (G::Graph)->
local V:= {$nops(op(3,G))}, A:= op(4,G) union~ `{}`~(V), s, v, U:= table([{}= {}]), Us;
    in V do
        for s in {indices}(U, 'nolist') do
            Us:= U[s]; U[s]:= evaln(U[s]);
            for v in V minus s do
                if (U[s union {v}]:= Us union A[v]) = V then return nops({s[], v}) fi
            od
        od
    od
:
 
> 

Bondage_number:=proc(g::Graph)
local dom,vn,edge_ofg,t,M,removeedge,H,i,hasNext, getNext,result;
dom:=Domination_number(g);
vn:=nops(op(3,g));
edge_ofg:=convert(Edges(g),list);
result:="Null";
for t from 1 to vn do
    M := Combination(vn, t);
    hasNext, getNext := ModuleIterator(M);
      while hasNext() do
        removeedge:= {seq(edge_ofg[i], i in getNext()+~1)};
        H := DeleteEdge(g, removeedge,inplace = false);
        if Domination_number(H)> dom then       
         result:=nops(removeedge);
         break t;
        end if;  
      end do;
end do:
return result;
end proc:
> 

ed:={{1, 2}, {1, 3}, {1, 4}, {1, 5}, {1, 7}, {1, 9}, {1, 11}, {2, 3}, {2, 5},
 {2, 6}, {2, 7}, {2, 12}, {2, 14}, {3, 4}, {3, 7}, {3, 8}, {3, 9}, {3, 16}, {4, 5},
 {4, 9}, {4, 10}, {4, 11}, {4, 17}, {5, 6}, {5, 11}, {5, 12}, {5, 19}, {6, 7}, {6, 12},{6, 13}, {6, 14}, {6, 21}, {7, 8}, {7, 14}, {7, 15}, {8, 9}, {8, 14}, {8, 15}, {8, 16}, {8, 22}, {9, 10}, {9, 16}, {9, 17}, {10, 11}, {10, 17}, {10, 18}, {10, 19}, {10, 23}, {11, 12},{11, 18}, {11, 19}, {12, 13}, {12, 19}, {12, 20}, {13, 14}, {13, 19}, {13, 20}, {13, 21}, {13, 24},{14, 15}, {14, 21}, {15, 16}, {15, 21}, {15, 22}, {15, 24}, {16, 17}, {16, 22}, {16, 23}, {17, 18},
{17, 22}, {17, 23}, {18, 19}, {18, 20}, {18, 23}, {18, 24}, {19, 20}, {20, 21}, {20, 23}, {20, 24},{21, 22}, {21, 24}, {22, 23}, {22, 24}, {23, 24}}:
g:=Graph(ed);

GRAPHLN(undirected, unweighted, [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24], Array(1..24, {(1) = {2, 3, 4, 5, 7, 9, 11}, (2) = {1, 3, 5, 6, 7, 12, 14}, (3) = {1, 2, 4, 7, 8, 9, 16}, (4) = {1, 3, 5, 9, 10, 11, 17}, (5) = {1, 2, 4, 6, 11, 12, 19}, (6) = {2, 5, 7, 12, 13, 14, 21}, (7) = {1, 2, 3, 6, 8, 14, 15}, (8) = {3, 7, 9, 14, 15, 16, 22}, (9) = {1, 3, 4, 8, 10, 16, 17}, (10) = {4, 9, 11, 17, 18, 19, 23}, (11) = {1, 4, 5, 10, 12, 18, 19}, (12) = {2, 5, 6, 11, 13, 19, 20}, (13) = {6, 12, 14, 19, 20, 21, 24}, (14) = {2, 6, 7, 8, 13, 15, 21}, (15) = {7, 8, 14, 16, 21, 22, 24}, (16) = {3, 8, 9, 15, 17, 22, 23}, (17) = {4, 9, 10, 16, 18, 22, 23}, (18) = {10, 11, 17, 19, 20, 23, 24}, (19) = {5, 10, 11, 12, 13, 18, 20}, (20) = {12, 13, 18, 19, 21, 23, 24}, (21) = {6, 13, 14, 15, 20, 22, 24}, (22) = {8, 15, 16, 17, 21, 23, 24}, (23) = {10, 16, 17, 18, 20, 22, 24}, (24) = {13, 15, 18, 20, 21, 22, 23}}), `GRAPHLN/table/1`, 0)

 (1)
> 

Domination_number(g)

4

 (2)
> 

Bondage_number(g);#Long time..
 

 

Bondage_number.mw

For arbitrarily given graph, it has been proved that determining its bondage number is NP-hard.

I always have a feeling that when removing subsets of edges, the symmetry of these edge subsets can be utilized. Above, we blindly traverse k-subsets of edges. If the graph is symmetric enough, perhaps we can reduce the amount of traversal. However, I have been struggling to understand how to express the symmetry of these edge subsets. Recently, I raised a similar question. 

 

 

 

error in f solve (AGM )...

Asked by: KIRAN SAJJAN 35
rootfinding fsolve numeric
June 02 2023
0  5

 

> 

restart; with(plots)
> 

PDEtools[declare]((F, T, G, H)(Y), prime = Y)

` F`(Y)*`will now be displayed as`*F

  

` T`(Y)*`will now be displayed as`*T

  

` G`(Y)*`will now be displayed as`*G

  

` H`(Y)*`will now be displayed as`*H

  

`derivatives with respect to`*Y*`of functions of one variable will now be displayed with '`

 (1)
> 

p1 := 0.1e-1; p2 := 0.3e-1; p3 := 0.5e-1; p := p1+p2+p3
> 

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

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

sigma2 := 0.210e-5; rs2 := 5700; ks2 := 25; cps2 := 523
> 

sigma3 := 6.30*10^7; rs3 := 10500; ks3 := 429; cps3 := 235
> 

sigma4 := 10^(-10); rs4 := 3970; ks4 := 40; cps4 := 765
> 

sigma5 := 1.69*10^7; rs5 := 7140; ks5 := 116; cps5 := 390
> 

sigma6 := 4.10*10^7; rs6 := 19300; ks6 := 318; cps6 := 129
> 

``
> 

M := 1; S1 := .5; A := 1; delta := 0.1e-2; g := .1; Gr := .5; betu := .5; bett := .5; Pr := 21; Ec := .5; bet := 1; S2 := .5; Rd := 1; Q := .1; Ra := .5; S := .1
> 

alp := .1
> 

``
> 

``
> 

B1 := 1+2.5*p+6.2*p^2; B2 := 1+13.5*p+904.4*p^2; B3 := 1+37.1*p+612.6*p^2; B4 := (ks1+2*kf-2*p*(kf-ks1))/(ks1+2*kf+p*(kf-ks1)); B5 := (ks2+3.9*kf-3.9*p*(kf-ks2))/(ks2+3.9*kf+p*(kf-ks2)); B6 := (ks3+4.7*kf-4.7*p*(kf-ks3))/(ks3+4.7*kf+p*(kf-ks3)); B7 := (ks4+2*kf-2*p*(kf-ks4))/(ks4+2*kf+p*(kf-ks4)); B8 := (ks5+3.9*kf-3.9*p*(kf-ks5))/(ks5+3.9*kf+p*(kf-ks5)); B9 := (ks6+4.7*kf-4.7*p*(kf-ks6))/(ks6+4.7*kf+p*(kf-ks6))
> 

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

a2 := 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
> 

``
> 

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
> 

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

a7 := 1-p1-p2-p3+p1*rs4/rf+p2*rs5/rf+p3*rs6/rf
> 

a8 := 1-p1-p2-p3+p1*rs4*cps4/(rf*cpf)+p2*rs5*cps5/(rf*cpf)+p3*rs6*cps6/(rf*cpf)
> 

a9 := B7*p1+B8*p2+B9*p3
> 

``
> 

a10 := 1+3*((p1*sigma4+p2*sigma5+p3*sigma6)/sigmaf-p1-p2-p3)/(2+(p1*sigma4+p2*sigma5+p3*sigma6)/((p1+p2+p3)*sigmaf)-((p1*sigma4+p2*sigma5+p3*sigma6)/sigmaf-p1-p2-p3))
> 

W := sum(b[i]*Y^i, i = 0 .. 3); Theta := sum(c[i]*Y^i, i = 0 .. 3); U := sum(d[i]*Y^i, i = 0 .. 2); Phi := sum(h[i]*Y^i, i = 0 .. 2)

Y^3*b[3]+Y^2*b[2]+Y*b[1]+b[0]

  

Y^3*c[3]+Y^2*c[2]+Y*c[1]+c[0]

  

Y^2*d[2]+Y*d[1]+d[0]

  

Y^2*h[2]+Y*h[1]+h[0]

 (2)
> 

F := a1*(1+1/bet)*(diff(W, `$`(Y, 2)))+a2*Ra*(diff(W, Y))+A-a5*M*W-S2*W^2+a2*Gr*Theta-S*betu*(W-U) = 0

9.1682928*Y*b[3]+3.0560976*b[2]+2.433571428*Y^2*b[3]+1.622380952*Y*b[2]+.8111904760*b[1]+1-1.346703274*b[3]*Y^3-1.346703274*b[2]*Y^2-1.346703274*b[1]*Y-1.346703274*b[0]-.5*(Y^3*b[3]+Y^2*b[2]+Y*b[1]+b[0])^2+.8111904760*c[3]*Y^3+.8111904760*c[2]*Y^2+.8111904760*c[1]*Y+.8111904760*c[0]+0.5e-1*d[2]*Y^2+0.5e-1*d[1]*Y+0.5e-1*d[0] = 0

 (3)
> 

T := (a4+Rd)*(diff(Theta, `$`(Y, 2)))+a3*Pr*Ra*(diff(Theta, Y))+Q*Theta+Pr*alp*S*bett*(Theta-Phi)+Pr*Ec*((1+1/bet)*a1*(diff(W, Y))^2+a5*M*W^2+(1+1/bet)*a1*S1*W^2+S2*W^3+S*betu*(W-U)) = 0

.205*c[1]*Y+.205*c[2]*Y^2+.205*c[3]*Y^3-.525*d[1]*Y-.525*d[2]*Y^2+.525*b[0]+.525*b[1]*Y+.525*b[2]*Y^2+.525*b[3]*Y^3+21.63764058*(Y^3*b[3]+Y^2*b[2]+Y*b[1]+b[0])^2-.105*h[0]+6.799682664*Y*c[3]+30.71903373*Y^2*c[3]+20.47935582*Y*c[2]-.105*Y^2*h[2]-.105*Y*h[1]-.525*d[0]+.205*c[0]+10.23967791*c[1]+2.266560888*c[2]+16.04451240*(3*Y^2*b[3]+2*Y*b[2]+b[1])^2+5.25*(Y^3*b[3]+Y^2*b[2]+Y*b[1]+b[0])^3 = 0

 (4)
> 

G := Ra*(diff(U, Y))+betu*(W-U) = 0

1.0*Y*d[2]+.5*d[1]+.5*b[3]*Y^3+.5*b[2]*Y^2-.5*d[2]*Y^2+.5*b[1]*Y-.5*d[1]*Y+.5*b[0]-.5*d[0] = 0

 (5)
> 

H := Ra*(diff(Phi, Y))+bett*(Theta-Phi) = 0

1.0*Y*h[2]+.5*h[1]+.5*c[3]*Y^3+.5*c[2]*Y^2-.5*Y^2*h[2]+.5*c[1]*Y-.5*Y*h[1]+.5*c[0]-.5*h[0] = 0

 (6)
> 

BCS := (D(W))(0) = 0, (D(Theta))(0) = 0, W(1) = -delta*(1+1/bet)*(D(W))(1), Theta(1) = 1+g*(D(Theta))(1), U(1) = -delta*(1+1/bet)*(D(W))(1), Phi(1) = 1+g*(D(Theta))(1)
> 

W := unapply(W(Y), Y)
> 

F := unapply(F(Y), Y)
> 

Theta := unapply(Theta(Y), Y)
> 

T := unapply(T(Y), Y)
> 

U := unapply(U(Y), Y)
> 

G := unapply(G(Y), Y)
> 

Phi := unapply(Phi(Y), Y)
> 

H := unapply(H(Y), Y)
> 

z1 := (D(W))(0) = 0

(D(W))(0) = 0

 (7)
> 

z2 := (D(Theta))(0) = 0

(D(Theta))(0) = 0

 (8)
> 

z3 := W(1) = -delta*(1+1/bet)*(D(W))(1)

b[3](1)+b[2](1)+b[1](1)+b[0](1) = -0.2e-2*(D(W))(1)

 (9)
> 

z4 := Theta(1) = 1+g*(D(Theta))(1)

c[3](1)+c[2](1)+c[1](1)+c[0](1) = 1+.1*(D(Theta))(1)

 (10)
> 

z5 := U(1) = -delta*(1+1/bet)*(D(W))(1)

d[2](1)+d[1](1)+d[0](1) = -0.2e-2*(D(W))(1)

 (11)
> 

z6 := Phi(1) = 1+g*(D(Theta))(1)

h[2](1)+h[1](1)+h[0](1) = 1+.1*(D(Theta))(1)

 (12)
> 

z7 := F(0)

1.+3.0560976*b[2](0)+.8111904760*b[1](0)-1.346703274*b[0](0)-.5*b[0](0)^2+.8111904760*c[0](0)+0.5e-1*d[0](0) = 0

 (13)
> 

z8 := T(0)

.525*b[0](0)+21.63764058*b[0](0)^2-.105*h[0](0)-.525*d[0](0)+.205*c[0](0)+10.23967791*c[1](0)+2.266560888*c[2](0)+16.04451240*b[1](0)^2+5.25*b[0](0)^3 = 0

 (14)
> 

z9 := G(0)

.5*d[1](0)+.5*b[0](0)-.5*d[0](0) = 0

 (15)
> 

z10 := H(0)

.5*h[1](0)+.5*c[0](0)-.5*h[0](0) = 0

 (16)
> 

z11 := F(1)

10.25516095*b[3](1)+3.331775278*b[2](1)-.5355127980*b[1](1)+1-1.346703274*b[0](1)-.5*(b[3](1)+b[2](1)+b[1](1)+b[0](1))^2+.8111904760*c[3](1)+.8111904760*c[2](1)+.8111904760*c[1](1)+.8111904760*c[0](1)+0.5e-1*d[2](1)+0.5e-1*d[1](1)+0.5e-1*d[0](1) = 0

 (17)
> 

z12 := T(1)

10.44467791*c[1](1)+22.95091671*c[2](1)+37.72371639*c[3](1)-.525*d[1](1)-.525*d[2](1)+.525*b[0](1)+.525*b[1](1)+.525*b[2](1)+.525*b[3](1)+21.63764058*(b[3](1)+b[2](1)+b[1](1)+b[0](1))^2-.105*h[0](1)-.105*h[2](1)-.105*h[1](1)-.525*d[0](1)+.205*c[0](1)+16.04451240*(3*b[3](1)+2*b[2](1)+b[1](1))^2+5.25*(b[3](1)+b[2](1)+b[1](1)+b[0](1))^3 = 0

 (18)
> 

z13 := G(1)

.5*d[2](1)+.5*b[3](1)+.5*b[2](1)+.5*b[1](1)+.5*b[0](1)-.5*d[0](1) = 0

 (19)
> 

z14 := H(1)

.5*h[2](1)+.5*c[3](1)+.5*c[2](1)+.5*c[1](1)+.5*c[0](1)-.5*h[0](1) = 0

 (20)
> 

NULL
> 

Z := fsolve([z1, z2, z3, z4, z5, z6, z7, z8, z9, z10, z11, z12, z13, z14], {b[0], b[1], b[2], b[3], c[0], c[1], c[2], c[3], d[0], d[1], d[2], h[0], h[1], h[2]})

(21)
> 

"F(Y):=eval(sum(b[i]*Y^i,i=0..5),Z); "

Error, invalid input: eval received NULL, which is not valid for its 2nd argument, eqns

  
> 

"T(Y):=eval(sum(c[i]*Y^i,i=0..5),Z); "

Error, invalid input: eval received NULL, which is not valid for its 2nd argument, eqns

  
> 

"G(Y):=eval(sum(d[i]*Y^i,i=0..5),Z); "

Error, invalid input: eval received NULL, which is not valid for its 2nd argument, eqns

  
> 

"H(Y):=eval(sum(h[i]*Y^i,i=0..5),Z); "

Error, invalid input: eval received NULL, which is not valid for its 2nd argument, eqns

  
> 

NULL
> 

plot(F(Y), Y = 0 .. 1, axes = boxed, color = green, thickness = 2, labels = [Y, W])

Error, (in plot) unexpected options: [9.1682928*Y(Y)*b[3](Y)+3.0560976*b[2](Y)+2.433571428*Y(Y)^2*b[3](Y)+1.622380952*Y(Y)*b[2](Y)+.8111904760*b[1](Y)+1-1.346703274*b[3](Y)*Y(Y)^3-1.346703274*b[2](Y)*Y(Y)^2-1.346703274*b[1](Y)*Y(Y)-1.346703274*b[0](Y)-.5*(b[3](Y)*Y(Y)^3+b[2](Y)*Y(Y)^2+b[1](Y)*Y(Y)+b[0](Y))^2+.8111904760*c[3](Y)*Y(Y)^3+.8111904760*c[2](Y)*Y(Y)^2+.8111904760*c[1](Y)*Y(Y)+.8111904760*c[0](Y)+0.5e-1*d[2](Y)*Y(Y)^2+0.5e-1*d[1](Y)*Y(Y)+0.5e-1*d[0](Y) = 0, Y = 0 .. 1]

  
> 

plot(T(Y), Y = 0 .. 1, axes = boxed, color = green, thickness = 2, labels = [Y, Theta])

Error, (in plot) unexpected options: [.205*c[1](Y)*Y(Y)+.205*c[2](Y)*Y(Y)^2+.205*c[3](Y)*Y(Y)^3-.525*d[1](Y)*Y(Y)-.525*d[2](Y)*Y(Y)^2+.525*b[0](Y)+.525*b[1](Y)*Y(Y)+.525*b[2](Y)*Y(Y)^2+.525*b[3](Y)*Y(Y)^3+21.63764058*(b[3](Y)*Y(Y)^3+b[2](Y)*Y(Y)^2+b[1](Y)*Y(Y)+b[0](Y))^2-.105*h[0](Y)+6.799682664*Y(Y)*c[3](Y)+30.71903373*Y(Y)^2*c[3](Y)+20.47935582*Y(Y)*c[2](Y)-.105*Y(Y)^2*h[2](Y)-.105*Y(Y)*h[1](Y)-.525*d[0](Y)+.205*c[0](Y)+10.23967791*c[1](Y)+2.266560888*c[2](Y)+16.04451240*(3*Y(Y)^2*b[3](Y)+2*Y(Y)*b[2](Y)+b[1](Y))^2+5.25*(b[3](Y)*Y(Y)^3+b[2](Y)*Y(Y)^2+b[1](Y)*Y(Y)+b[0](Y))^3 = 0, Y = 0 .. 1]

  
> 

plot(G(Y), Y = 0 .. 1, axes = boxed, color = green, thickness = 2, labels = [Y, W])

Error, (in plot) unexpected options: [1.0*Y(Y)*d[2](Y)+.5*d[1](Y)+.5*b[3](Y)*Y(Y)^3+.5*b[2](Y)*Y(Y)^2-.5*d[2](Y)*Y(Y)^2+.5*b[1](Y)*Y(Y)-.5*d[1](Y)*Y(Y)+.5*b[0](Y)-.5*d[0](Y) = 0, Y = 0 .. 1]

  
> 

plot(H(Y), Y = 0 .. 1, axes = boxed, color = green, thickness = 2, labels = [Y, Theta])

Error, (in plot) unexpected options: [1.0*Y(Y)*h[2](Y)+.5*h[1](Y)+.5*c[3](Y)*Y(Y)^3+.5*c[2](Y)*Y(Y)^2-.5*Y(Y)^2*h[2](Y)+.5*c[1](Y)*Y(Y)-.5*Y(Y)*h[1](Y)+.5*c[0](Y)-.5*h[0](Y) = 0, Y = 0 .. 1]

  

NULL


 

Download AGM.mw
 

 

Error, (in dsolve/numeric/process_input) input sys...

Asked by: Saha 70
ode syntax homework
June 01 2023
0  3

I have solved the  system with dsolve (using numeric option).But i got this Error.

Here is my code.

RK_EC.mw

does mserver.exe use threads? How to turn that off...

Asked by: nm 8552
multithreading mserver
June 01 2023
0  1

In trying to find why my Maple program encouter long delays when using timelimit and hangs, I found from windows 10 task manager that it uses thread and the thread can hang on network IO. I have no idea why mserver.exe is doing network IO for in the middle of timelimit(). 

But I think this has something to do with the problems I am seeing (reference).

The first thing I noticed is the network panel has check for update ON. So I turned that off.

I do not know if this was causing the problem, where Maple in middle of computation was trying to check for an update over the network or not.

But my question here is this: Does mserver.exe uses threads under the cover? If so, is there an option to turn this off? I.e tell Maple server.exe not to use threads at all?

I was to see if this is the cause or not. I see on the help on kernel options it says

And when I do   kernelopts(multithreaded)  it says  true

But how to turn this off? When I type

kernelopts(multithreaded=false)
Error, kernelopts cannot set multithreaded value

Is it enough to tell it to use ONE thread only then? Like this

kernelopts(gcmaxthreads=1)
                            numcpus

kernelopts(gcmaxthreads)
1

If one is not able turn multhreading off, will setting gcmaxthreads=1 have same effect or is there a better way to do these things?

My code does not do any mutlithreading. So I do not need it.

Maple 2023 on windows 10

Is it possible to permanently turn off threading? I looked at options and see no such option.

 

 

 

How can I calculate the domination number of a gra...

Asked by: lcz 994
May 31 2023
2  10

In graph theory, a dominating set for a graph G is a subset D of its vertices, such that any vertex of G is either in D, or has a neighbor in D. The domination number γ(G) is the number of vertices in a smallest dominating set for G. The domination number is a well-known parameter in graph theory. But I am unable to find a built-in function in Maple to calculate the domination number of a graph. Did I miss something?

For example, I would like to calculate the domination number of the following graph.

ed:={{1, 2}, {1, 3}, {1, 4}, {1, 5}, {1, 7}, {1, 9}, 
{1, 11}, {2, 3}, {2, 5}, {2, 6}, {2, 7}, {2, 12}, {2, 14},
 {3, 4}, {3, 7}, {3, 8}, {3, 9}, {3, 16}, {4, 5}, {4, 9}, 
{4, 10}, {4, 11}, {4, 17}, {5, 6}, {5, 11}, {5, 12}, {5, 19}, 
{6, 7}, {6, 12}, {6, 13}, {6, 14}, {6, 21}, {7, 8}, {7, 14}, 
{7, 15}, {8, 9}, {8, 14}, {8, 15}, {8, 16}, {8, 22}, {9, 10},
 {9, 16}, {9, 17}, {10, 11}, {10, 17}, {10, 18}, {10, 19}, 
{10, 23}, {11, 12}, {11, 18}, {11, 19}, {12, 13}, {12, 19}, 
{12, 20}, {13, 14}, {13, 19}, {13, 20}, {13, 21}, {13, 24}, 
{14, 15}, {14, 21}, {15, 16}, {15, 21}, {15, 22}, {15, 24}, 
{16, 17}, {16, 22}, {16, 23}, {17, 18}, {17, 22}, {17, 23}, 
{18, 19}, {18, 20}, {18, 23}, {18, 24}, {19, 20}, {20, 21}, 
{20, 23}, {20, 24}, {21, 22}, {21, 24}, {22, 23}, {22, 24}, {23, 24}};
g:=Graph(ed);

 

 

Inconsistent behavior with dsolve/numeric and NLPS...

Asked by: Dr. Venkat Subramanian 386
optimization ode dsolve
May 31 2023
0  34

dsolve/numeric + NLPSolve shows inconsistent behavior. This combination is important for parameter estimation and optimal control. Can anyone fix this? Hopefully, I am not making a mistake.

> 

restart:
> 

 

Test code written by Dr. Venkat Subramanian at UT Austin, 05/31/2023. This code uses CVP approach (piecwise constant) to perform optimal control. NLPSolve combination with dsolve numeric parametric form is buggy and fails for some values of nvar, and works for some values of nvar. Ideally increasing nvar should show convergence with respect to the objective function.

> 

restart:
> 

Digits:=15;

15

> 

 

> 

eqodes:=[diff(ca(t),t)=-(u+u^2/2)*1.0*ca(t),diff(cb(t),t)=1.0*u*ca(t)-0.1*cb(t)];

[diff(ca(t), t) = -1.0*(u+(1/2)*u^2)*ca(t), diff(cb(t), t) = 1.0*u*ca(t)-.1*cb(t)]

> 

soln:=dsolve({op(eqodes),ca(0)=alpha,cb(0)=beta},type=numeric,'parameters'=[alpha,beta,u],compile=true,savebinary=true):
> 

 

> 

 

> 

ss:=proc(x)
interface(warnlevel=0):
#if  type(x[1],numeric)
if  type(x,Vector)
then

local z1,n1,i,c10,c20,dt,u;
global soln,nvar;
dt:=evalf(1.0/nvar):
c10:=1.0:c20:=0.0:
for i from 1 to nvar do
u:=x[i]:
soln('parameters'=[c10,c20,u]):
z1:=soln(dt):
c10:=subs(z1,ca(t)):c20:=subs(z1,cb(t)):
od:
-c20;
 else 'procname'(args):

end if:

end proc:
> 

 

> 

nvar:=3;#code works for nvar:=3, but not for nvar:=2

3

> 

ic0:=Vector(nvar,[seq(0.1,i=1..nvar)],datatype=float[8]);

Vector(3, {(1) = .1000000000000000, (2) = .1000000000000000, (3) = .1000000000000000})

> 

 

> 

ss(ic0);

HFloat(-0.09025793628011441)

> 

bl := Vector(nvar,[seq(0.,i=1..nvar)]);bu := Vector(nvar,[seq(5.,i=1..nvar)]);

Vector(3, {(1) = 0., (2) = 0., (3) = 0.})

Vector[column](%id = 36893491136404053036)

> 

infolevel[all]:=1;

1

> 

Optimization:-NLPSolve(nvar,ss,[],initialpoint=ic0,[bl,bu]);

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 3

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

[-.531453381886523, Vector(3, {(1) = .8114345197312305, (2) = 1.189413022736622, (3) = 2.427689509710160})]

> 

nvar:=2;#code works for nvar:=3, but not for nvar:=2

2

> 

ic0:=Vector(nvar,[seq(0.1,i=1..nvar)],datatype=float[8]);

Vector(2, {(1) = .1000000000000000, (2) = .1000000000000000})

> 

 

> 

ss(ic0);

HFloat(-0.09025793011810783)

> 

bl := Vector(nvar,[seq(0.,i=1..nvar)]);bu := Vector(nvar,[seq(5.,i=1..nvar)]);

Vector(2, {(1) = 0., (2) = 0.})

Vector[column](%id = 36893491136437818540)

> 

infolevel[all]:=1;

1

> 

Optimization:-NLPSolve(nvar,ss,[],initialpoint=ic0,[bl,bu]);

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

Error, (in Optimization:-NLPSolve) no improved point could be found

> 

 

> 

nvar:=5;#code works for nvar:=3, but not for nvar:=2

5

> 

ic0:=Vector(nvar,[seq(0.1,i=1..nvar)],datatype=float[8]);

Vector(5, {(1) = .1000000000000000, (2) = .1000000000000000, (3) = .1000000000000000, (4) = .1000000000000000, (5) = .1000000000000000})

> 

 

> 

ss(ic0);

HFloat(-0.09025792639212991)

> 

bl := Vector(nvar,[seq(0.,i=1..nvar)]);bu := Vector(nvar,[seq(5.,i=1..nvar)]);

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

Vector[column](%id = 36893491136472713804)

> 

infolevel[all]:=1;

1

> 

Optimization:-NLPSolve(nvar,ss,[],initialpoint=ic0,[bl,bu]);

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 5

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

[-.537338871783244, Vector(5, {(1) = .7435767224059456, (2) = .9013440906589676, (3) = 1.148772394841535, (4) = 1.629079877401040, (5) = 3.222801229872320})]

> 

 


 

Download test.mw

First 74 75 76 77 78 79 80 Last Page 76 of 2308