janhardo

700 Reputation

12 Badges

11 years, 41 days

MaplePrimes Activity


These are replies submitted by janhardo

LAD_complex := proc(pde, Uxt::function(name), inx, t::name, N::posint)
    option `Complex LAD Procedure v0.1`;
    local Ru, u, pde1, Ut, q1, q2, scale, term, T, RRu, NNu, A, s, lambda, vars, nvars, vvars, NNv, i, v, vals, xt, n, j, Ai, soln, z;

    # Step 1: Prevent time variable in the initial condition
    if has(inx, t) then error "initial condition cannot contain %1", t end if;

    # Step 2: Normalize the PDE into standard form
    if type(pde, `=`) then
        pde1 := lhs(pde) - rhs(pde)
    else
        pde1 := expand(pde)
    end if;
    if not type(pde1, `+`) then error "PDE must be a sum of terms" end if;

    # Step 3: Identify the time derivative term
    Ut := diff(Uxt, t);
    q1, q2 := selectremove(has, pde1, Ut);
    if q1 = 0 then error "time derivative of %1 not found", Uxt end if;

    # Step 4: Verify time term is linear in I * dU/dt
    scale := coeff(q1, Ut);
    if indets(scale) minus {I} <> {} then error "time term must be linear in I * dU/dt" end if;

    # Step 5: Isolate the right-hand side and normalize it
    q2 := expand(-q2 / scale);

    # Step 6: Separate linear and nonlinear parts
    RRu, NNu := selectremove(term -> not has(eval(term, Uxt = T*Uxt)/T, T), q2 + z);
    NNu := NNu - z;

    # Step 7: Initialize series terms and variables
    xt := op(Uxt);
    u[0] := inx;
    soln := inx;

    vars := select(has, [op(indets(NNu))], Uxt);
    nvars := nops(vars);
    vvars := [seq(v[i], i = 1 .. nvars)];
    NNv := subs(vars =~ vvars, NNu);

    vals[0] := simplify(eval(subs(Uxt = u[0], vars))) assuming real;
    Ai := eval(subs(vvars =~ vals[0], NNv));
    A[0] := simplify(Ai) assuming real;
    Ru := simplify(eval(subs(Uxt = u[0], RRu)));

    # Step 8: Iterative expansion in time
    for n to N do
        u[n] := inttrans:-invlaplace(inttrans:-laplace(Ru + A[n - 1], t, s)/s, s, t);
        soln := soln + u[n];
        vals[n] := simplify(eval(subs(Uxt = u[n], vars)));
        A[n] := simplify(eval(diff(subs({seq(v[j] = add(vals[i][j]*lambda^i, i = 0 .. n), j = 1 .. nvars)}, NNv), lambda $ n)/n!, lambda = 0));
        Ru := simplify(eval(subs(Uxt = u[n], RRu)));
    end do;

    simplify(soln) assuming real;
end proc:
 

pde := diff(U(x, t), t) + U(x, t)^2*diff(U(x, t), x);
infolevel[LAD] := 2;
inx := 3*x;
N := 10;
LAD_complex(pde, U(x, t), inx, t, 10);

diff(U(x, t), t)+U(x, t)^2*(diff(U(x, t), x))

 

2

 

3*x

 

10

 

175692092397588*t^10*x^11-5650915252554*t^9*x^10+184670433090*t^8*x^9-6155681103*t^7*x^8+210450636*t^6*x^7-7440174*t^5*x^6+275562*t^4*x^5-10935*t^3*x^4+486*t^2*x^3-27*t*x^2+3*x

(1)

pde1 := diff(U(x, t), t)*I + a[1]*diff(U(x, t), x, x, x)*I + b*U(x, t)*conjugate(U(x, t))*U(x, t);
init := 3*x;
LAD_complex(pde1, U(x, t), init, t, 4);

(diff(U(x, t), t))*I+I*a[1]*(diff(diff(diff(U(x, t), x), x), x))+b*U(x, t)^2*conjugate(U(x, t))

 

3*x

 

(1/8)*(-708588*(((5/81)*I)*b*t*abs(x)^14+(-(1/1458)*t*a[1]+2/729)*abs(x)^12-(80/243)*t*a[1]*abs(x)^10-(20/81)*t*a[1]*abs(x)^8+(I*b*t*x^4-(4/243)*a[1]*(x^2+20)*t+(8/243)*x^2)*x^9)*abs(x)*b^2*t^3*a[1]*signum(1, x)+271188*((3/124)*b^4*t^4*x^9+I*t^3*(-(1/93)*x^3+t*a[1])*x^4*b^3+(125/558)*t^2*(a[1]^2*(x^2+(36/5)*x+54/5)*t^2+(44/125)*t*x*a[1]-(2/125)*x^4)*x*b^2-((2/837)*I)*(-(1/3)*x^3+t*a[1])*t*b+(2/22599)*x)*x^9)/x^9

(2)

Download Uitbreiding_LAD_procedure_naar_complex_12-6-2025mprimes.mw

@dharr 
What you could do is approach it more systematically and use the LAD procedure as a basis.
Extend the PDE:

by adding the simplest new term or expression.
A possible check in the form of: 'Extrapolate to an infinite number of terms' probably won't exist for that extended PDE.




Yes, you can choose initial conditions arbitrarily — but not randomly if your goal is to analyze specific behavior.

But Beware:

  • Choosing ICs too far from equilibrium can lead to misleading visuals, especially if you're interested in local behavior near a critical point.

  • In linear systems, every IC leads to a unique trajectory because the system is deterministic and satisfies the existence and uniqueness theorem.

In this example, it are 4 points of a square around the critical point.

Adding more points and lines,circle 


 

protected names O and I   , that's why using O_ , I_   

with(plots); alpha := (1/4)*Pi; a := 2^alpha*exp(I*alpha); b := 2^(alpha+1)*exp(I*(alpha+(1/2)*Pi)); h := 1/(1/a+1/b); O_ := 0; A := a; B := b; H := h; I_ := (O_+H)*(1/2); J := (H+B)*(1/2); toCoords := proc (z) options operator, arrow; [Re(evalf(z)), Im(evalf(z))] end proc; pO := toCoords(O_); pA := toCoords(A); pB := toCoords(B); pH := toCoords(H); pI := toCoords(I_); pJ := toCoords(J); triangle := plot([pO, pA, pB, pO], color = blue, thickness = 2); points := pointplot([pO, pA, pB, pH, pI, pJ], symbol = solidcircle, color = red, symbolsize = 15); xaxis := textplot([[2, -.3, "Re(z)"]], font = [TIMES, BOLD, 12]); yaxis := textplot([[-.3, 2, "Im(z)"]], font = [TIMES, BOLD, 12]); makeLabel := proc (name, coords) textplot([coords[1]+0.5e-1, coords[2], cat(name, "(", sprintf("%.2f", coords[1]), ",", sprintf("%.2f", coords[2]), ")")]) end proc; labelO := makeLabel("O", pO); labelA := makeLabel("A", pA); labelB := makeLabel("B", pB); labelH := makeLabel("H", pH); labelI := makeLabel("I", pI); labelJ := makeLabel("J", pJ); display([triangle, points, xaxis, yaxis, labelO, labelA, labelB, labelH, labelI, labelJ], scaling = constrained, axes = normal, title = "Extended triangle OAB with points H, I, J", labels = ["Re(z)", "Im(z)"])

 

NULL


 

Download complexe_plotdriehoek_en_cirkel10-6-2025mprimes.mw

@jalal 
Hi,
a example to start with 

A := 1 + I;
B := 4 + 2*I;
C := 2 + 5*I;
punt := z -> [Re(z), Im(z)];
driehoek := plot([punt(A), punt(B), punt(C), punt(A)], color = blue, thickness = 2);
punten := pointplot([punt(A), punt(B), punt(C)], symbol = solidcircle, color = red, symbolsize = 15);
display([driehoek, punten], scaling = constrained, title = "Triangle in complex plane", axes = boxed);

It seems that the array notation is not suited for pattern matching in Maple ?

@jalal 
More intuiative with complex numbers ...

corrected code

restart;
with(DEtools);
with(plots);
with(LinearAlgebra);
f := y -> y^2*(1 - y^2);
sys := [diff(x(t), t) = -x(t), diff(y(t), t) = f(y(t))];
trajectories := [[x(0) = 1, y(0) = 1.5], [x(0) = -1, y(0) = 0.5], [x(0) = 0.5, y(0) = -1.5]];
phase_plot := DEplot(sys, [x(t), y(t)], t = 0 .. 10, trajectories, x = -2 .. 2, y = -2 .. 2, arrows = medium, dirfield = [20, 20], linecolor = blue, title = "Phase Portrait with Equilibrium Points");
equilibrium_points := pointplot([[0, 0], [0, 1], [0, -1]], symbol = solidcircle, color = black, symbolsize = 15);
display([phase_plot, equilibrium_points], scaling = constrained);
print("\n### DETAILED STABILITY ANALYSIS ###");
F := [-x, -y^4 + y^2];
J := Matrix([[diff(F[1], x), diff(F[1], y)], [diff(F[2], x), diff(F[2], y)]]);
print("\n1. Analysis for (0,0):");
J0 := eval(J, {x = 0, y = 0});
print("Jacobian:", J0);
print("Eigenvalues:", Eigenvalues(J0));
print("Conclusion:");
print("   - Eigenvalues: -1 and 0");
print("   - Non-hyperbolic point");
print("   - dy/dt &approx; y^2 > 0 for y &approx; 0");
print("   - &Implies; (0,0) is UNSTABLE");
print("\n2. Analysis for (0,1):");
J1 := eval(J, {x = 0, y = 1});
print("Jacobian:", J1);
print("Eigenvalues:", Eigenvalues(J1));
print("Conclusion:");
print("   - Eigenvalues: -1 and -2");
print("   - Both negative &Implies; STABLE NODE");
print("   - All nearby trajectories converge to this point");
print("\n3. Analysis for (0,-1):");
Jm1 := eval(J, {x = 0, y = -1});
print("Jacobian:", Jm1);
print("Eigenvalues:", Eigenvalues(Jm1));
print("Conclusion:");
print("   - Eigenvalues: -1 and +2");
print("   - Mixed eigenvalues &Implies; SADDLE POINT");
print("   - Stable in x-direction, unstable in y-direction");
print("   - &Implies; (0,-1) is UNSTABLE");
print("\n### VISUAL INTERPRETATION ###");
print("1. (0,0): Arrows point away in y-direction - unstable");
print("2. (0,1): All arrows point toward this point - stable");
print("3. (0,-1): Arrows approach in x-direction but diverge in y-direction - saddle point");
print("\n### SYSTEM BEHAVIOR SUMMARY ###");
print("The system exhibits:");
print("- One stable equilibrium at (0,1)");
print("- One unstable equilibrium at (0,0)");
print("- One saddle point at (0,-1)");
print("- Bistability between y=1 and y=-1 states");
print("- All x-values decay to 0 over time");
print("- y-direction determines long-term behavior");

Never trust ai , because node ( 0,-1) is not stable , when i examine it 

@salim-barzani 
You need help from a specialist , i don't know what you are after for.? 

@salim-barzani 

 

restart;
expr := abs(U[i](x, t))^n*diff(U[i](x, t), x);
expr_rewritten := subs(abs(U[i](x, t))^n = (U[i](x, t)*conjugate(U[i](x, t)))^(n/2), expr);
print(expr_rewritten);

"in here i want to seperate abs(U[i](x,t))^n=U[i](x,t)^n/2*conjugate(U[i](x,t))"   
I am guessing here, no direction where you are after for?

 


=======================================================================






 

 

 

 

1 2 3 4 5 6 7 Last Page 3 of 73