restart;

A := -x*(x - 4*exp(x/2) + 2);
B := x*sqrt((-8*x - 16)*exp(x/2) + x^2 + 4*x + 16*exp(x) + 4);

simplify(A-B) assuming x::real;
 

A2 := expand(A^2);
B2 := expand(B^2);

A2 := expand(A^2)assuming x::real;
B2 := expand(B^2)assuming x::real;

simplify(A2 - B2);

is(simplify(A2) = simplify(B2)) ;

solve ( x^4 - 8*x^3*exp(x/2) + 4*x^3 + 16*x^2*exp(x/2)^2 - 16*x^2*exp(x/2) + 4*x^2, x) assuming x::real;

solve( -8*x^3*exp(x/2) - 16*x^2*exp(x/2) + x^4 + 4*x^3 + 16*x^2*exp(x) + 4*x^2,x)assuming x::real;

restart;
with(PDEtools):
with(LinearAlgebra):
with(SolveTools):
undeclare(prime):
alias(F=F(x,y,z,t), G=G(x,y,z,t)):

ND := proc(F, G, U)
  local v, w, f, g, a, i, newv;
  v := [op(F)];
  newv := [];
  for i from 1 to nops(v) do
    if v[i] in U then
      newv := [op(newv), -v[i]];
    else
      newv := [op(newv), v[i]];
    end if;
  end do;
  f := op(0, F);
  g := op(0, G);
  a := diff(f(op(newv))*g(op(v)), U);
  convert(subs(newv=~v, a), diff);
end proc:

 

ND(sin(x)*cos(y)*exp(z-t), cos(x)*sin(y)*exp(z+t), [x,y,z,t]);

ND := proc(F, G, U)
  local v, w, f, g, a, i, newv, result;
  
  v := [op(F)];     # Haal variabelenlijst uit F
  newv := [];       # Maak lege lijst
  
  for i from 1 to nops(v) do
    if v[i] in U then
      newv := [op(newv), -v[i]];  # Spiegelen
    else
      newv := [op(newv), v[i]];   # Anders normaal houden
    end if;
  end do;
  
  f := op(0, F);    # Functienaam F
  g := op(0, G);    # Functienaam G
  
  a := diff(f(op(newv)) * g(op(v)), U);  # Differentieer aangepast product
  result := convert(subs(newv=~v, a), diff);  # Zet gespiegelde terug
  
  return simplify(result);  # >>> AUTOMATISCH VEREENVOUDIGEN! <<<
end proc:

ND(sin(x)*cos(y)*exp(z-t), cos(x)*sin(y)*exp(z+t), [x,y,z,t]);

ND(F, G, [x]);
ND(F, G, [t]);

ND(F(x,y,z,t), G(x,y,z,t), [x]);
ND(F(x,y,z,t), G(x,y,z,t), [y]);
ND(F(x,y,z,t), G(x,y,z,t), [z]);
ND(F(x,y,z,t), G(x,y,z,t), [t]);
ND(F(x,y,z,t), G(x,y,z,t), [x,y]);
ND(F(x,y,z,t), G(x,y,z,t), [x,y,z]);
ND(F(x,y,z,t), G(x,y,z,t), [x$2, y, z, t]);

ND(F(x,y,t), G(x,y,t), [x]);
ND(F(x,y,z,t), G(x,y,z,t), [y]);
ND(F(x,y,z,t), G(x,y,z,t), [z]);
ND(F(x,y,z,t), G(x,y,z,t), [t]);
ND(F(x,y,z,t), G(x,y,z,t), [x,y]);
ND(F(x,y,z,t), G(x,y,z,t), [x,y,z]);
ND(F(x,y,t), G(x,y,t), [x$2, y, t]);

restart;

# Definities
g := a1*x + a2*y + a3*z + a4*t + a5:
h := a6*x + a7*y + a8*z + a9*t + a10:
f := g^2 + h^2 + a11:

# Afgeleiden
fx := diff(f, x): fy := diff(f, y): fz := diff(f, z): ft := diff(f, t):
fxx := diff(fx, x): fxy := diff(fx, y): fxz := diff(fx, z): fyy := diff(fy, y):
fxxx := diff(fxx, x): fxxy := diff(fxx, y): fxxxy := diff(fxxx, y):
fxxxx := diff(fxxx, x): fxxxxx := diff(fxxxx, x): fxxxxxx := diff(fxxxxx, x):

# PDE
PDE := 9*(f*diff(f, x, t) - fx*ft)
    -5*(f*diff(f, x, x, x, y) - 3*fx*diff(f, x, x, y) + 3*fxx*fxy - fxxx*fy)
    + (f*diff(f, x, x, x, x, x, x) - 6*fx*diff(f, x, x, x, x, x) + 15*fxx*diff(f, x, x, x, x) - 10*(fxxx)^2)
    -5*(f*fyy - fy^2)
    +alpha*(f*fxx - fx^2)
    +beta*(f*fxy - fx*fy)
    +gamma*(f*fxz - fx*fz):

PDE := expand(PDE): PDE := simplify(PDE):

# Coëfficiënten pakken
coeff1 := coeff(PDE, x, 1):
coeff2 := coeff(PDE, y, 1):
coeff3 := coeff(PDE, z, 1):
coeff4 := coeff(PDE, t, 1):
coeff_const := eval(PDE, {x=0, y=0, z=0, t=0}):

# Stelsel van vergelijkingen
Eqns := {coeff1=0, coeff2=0, coeff3=0, coeff4=0, coeff_const=0}:
# Stelsel
##Eqns := {coeff_x=0, coeff_y=0, coeff_z=0, coeff_t=0, coeff_const=0}:

# Oplossen
Sol := solve(Eqns, {a4, a9, a11}):
Sol;

restart:

DrawOrderedLines := proc(shape1::list, isolated1::list, shape2::list, isolated2::list)
  local i, p1, p2, allplots, colorlist, polyline, isolateplot, j:
  uses plots, plottools;
  colorlist := [blue, red]:
  allplots := []:

  for i from 1 to 2 do
    if i = 1 then
      polyline := shape1:
      isolateplot := isolated1:
    else
      polyline := shape2:
      isolateplot := isolated2:
    end if:

    # Verbind punten in exacte volgorde
    if nops(polyline) >= 2 then
      for j from 1 to nops(polyline) - 1 do
        p1 := polyline[j]:
        p2 := polyline[j+1]:
        allplots := [op(allplots), line(p1, p2, color=colorlist[i], thickness=2)]:
      end do:
    end if:

    # Teken losse geïsoleerde punten
    if nops(isolateplot) > 0 then
      allplots := [op(allplots), pointplot(isolateplot, color=colorlist[i], symbol=solidcircle, symbolsize=15)]:
    end if:
  end do:

  return display(allplots, scaling=constrained, title="");
end proc:

fig1 := [[0, 0], [4, 0], [4, 3], [0, 3], [5, 5],[8,1],[8,8],[-5,8]]:
isolated1 := [[0,3],[2, 1],[-5,8]]:
fig2 := [[6, 6], [8, 6], [8, 8], [6, 8], [9, 9]]:
isolated2 := [[7, 7]]:

DrawOrderedLines(fig1, isolated1, fig2, isolated2);

restart:
with(DEtools):

ExactODEGroupingInteractive := proc(ODE, x, y)
  local lhs_expr, M, N, M_orig, N_orig, dM_dy, dN_dx, matched,
        intM, u_partial, du_dy, gprime, g, u_final, terms,
        dx_term, dy_term, ODE_classical, advice;

  print("�� STEP 0: Start by reading the differential equation in differential form:");
  lhs_expr := lhs(ODE);
  print("   → Equation provided: ", ODE);

  print("\n�� STEP 0B: Try to classify the ODE using Maple’s odeadvisor");

  # Extract M and N
  terms := [op(lhs_expr)];
  dx_term := select(t -> has(t, Dx), terms)[1];
  dy_term := select(t -> has(t, Dy), terms)[1];
  M := simplify(dx_term / Dx);
  N := simplify(dy_term / Dy);

  # Save for later
  M_orig := M:
  N_orig := N:

  # Substitute y = y(x)
  M := subs(y = y(x), M):
  N := subs(y = y(x), N):

  print("   We now rewrite the ODE in the form:");
  print("     M(x,y)*dy/dx + N(x,y) = 0");
  ODE_classical := M + N*diff(y(x), x) = 0;
  print("   → Rewritten ODE: ", ODE_classical);

  print("   Ask Maple’s odeadvisor for classification:");
  advice := odeadvisor(ODE_classical, y(x));
  print("�� Maple's odeadvisor says: ", advice);

  if member(_exact, advice) then
    print("✅ Great! Maple recognizes this ODE as exact.");
  else
    print("⚠️ Maple did NOT classify this as exact.");
    print("   → We will verify exactness ourselves using partial derivatives.");
  end if;

  print("\n�� STEP 1: Extract the functions M(x, y) and N(x, y) from the equation:");
  print("   → M(x, y) = ", M_orig);
  print("   → N(x, y) = ", N_orig);

  print("\n�� STEP 2: Check if the ODE is exact by comparing partial derivatives:");
  print("   → Compute ∂M/∂y and ∂N/∂x:");
  dM_dy := diff(M_orig, y):
  dN_dx := diff(N_orig, x):
  print("     ∂M/∂y = ", dM_dy);
  print("     ∂N/∂x = ", dN_dx);
  matched := simplify(dM_dy - dN_dx);
  if matched <> 0 then
    return "❌ This equation is NOT exact. Grouping method is not applicable.";
  else
    print("✅ Yes! The equation is exact.");
  end if;

  print("\n�� STEP 3: Integrate M(x, y) with respect to x (treat y as constant):");
  print("   → ∫ M(x, y) dx = ");
  intM := int(M_orig, x):
  print("     = ", intM, " + g(y)");

  print("\n�� STEP 4: Differentiate that result with respect to y:");
  print("   → ∂/∂y of the integral:");
  du_dy := diff(intM, y):
  print("     ∂(partial result)/∂y = ", du_dy);

  print("\n�� STEP 5: Now solve for g'(y):");
  print("   → g'(y) = N(x, y) - ∂/∂y(previous result)");
  gprime := simplify(N_orig - du_dy):
  print("     g'(y) = ", gprime);

  print("\n�� STEP 6: Integrate g'(y) with respect to y:");
  print("   → ∫ g'(y) dy = ");
  g := int(gprime, y):
  print("     g(y) = ", g);

  print("\n�� STEP 7: Add both parts to find the general solution u(x, y):");
  u_final := simplify(intM + g):
  print("     u(x, y) = ", u_final);

  print("\n�� FINAL RESULT: The general solution to the ODE is:");
  return u_final = 'C';

end proc:

# Example call
ExactODEGroupingInteractive(
  (6*x*y^2 + 4*x^3*y)*Dx + (6*x^2*y + x^4 + exp(y))*Dy = 0,
  x, y);

ExactODEGroupingInteractive(
  (3*x*y^2 + 4*x^3*y)*Dx + (6*x^2*y + x^4 + exp(y))*Dy = 0,
  x, y);

complex_map_plot := proc(f, zcurveLijst, t_range,
                         {kleurenset := ["blue","red","green","magenta","orange","cyan"],
                          lijnstijl := solid, resolutie := 200,
                          domein_naam := "Krommes in z-vlak", toon_label := false})
    local i, t, zcurve, ReZ, ImZ, ReW, ImW, zplots, wplots, kleur, zplot, wplot, annotatie, labeltekst;
    uses plots, plottools, InertForm;

    zplots := [];  wplots := [];

    for i to nops(zcurveLijst) do
        zcurve := zcurveLijst[i];
        kleur := kleurenset[i mod nops(kleurenset) + 1];

        ReZ := t -> Re(zcurve(t));
        ImZ := t -> Im(zcurve(t));
        ReW := t -> Re(f(zcurve(t)));
        ImW := t -> Im(f(zcurve(t)));

        zplot := plot([ReZ(t), ImZ(t), t = t_range],
                      color = kleur, linestyle = lijnstijl, thickness = 2,
                      scaling = constrained, axes = boxed, numpoints = resolutie,
                      labels = ["Re(z)", "Im(z)"]);

        if toon_label then
            labeltekst := typeset("z(t) = ", InertForm:-Display(zcurve(t)));
            annotatie := textplot([ReZ((rhs(t_range) + lhs(t_range))/2),
                                   ImZ((rhs(t_range) + lhs(t_range))/2),
                                   labeltekst], font = [TIMES, BOLD, 12]);
            zplot := display(zplot, annotatie);
        end if;

        wplot := plot([ReW(t), ImW(t), t = t_range],
                      color = kleur, linestyle = lijnstijl, thickness = 2,
                      scaling = constrained, axes = boxed, numpoints = resolutie,
                      labels = ["Re(w)", "Im(w)"]);

        zplots := [op(zplots), zplot];
        wplots := [op(wplots), wplot];
    end do;

    display(Array([display(zplots, title = cat("z-vlak: ", domein_naam)),
                   display(wplots, title = "w-vlak: Beeld onder f(z)")]), insequence = false);
end proc:

f := z -> exp(z);
zcurves := [t -> -1 + I*t, t -> 0 + I*t, t -> 1 + I*t, t -> 2 + I*t];
complex_map_plot(f, zcurves, -Pi..Pi, toon_label = false, domein_naam = "Verticale lijnen Re(z) = -1,0,1,2");

f := z -> exp(z);
z_driehoek := proc(t) piecewise(
  t <= 1, 2*t,
  t <= 2, 2 + (t-1)*(-1+I),
  t <= 3, 1 + I + (t-2)*(-1 - I)
) end proc;

complex_map_plot(f, [z_driehoek], 0..3, toon_label=true, domein_naam="Driehoek");

zcurves := [
    t -> t - 2*I,
    t -> t - I,
    t -> t,
    t -> t + I,
    t -> t + 2*I
];

f := z -> exp(z);
complex_map_plot(f, zcurves, -2..2, domein_naam = "Horizontale lijnen Im(z) = -2..2", toon_label = false);