I want to differentiate NAA with respect to U(L,J+1), how to do it.

2 visualizzazioni (ultimi 30 giorni)
%---------------------SPLINE IN COMPRESSION------------------------------------------------
% Numerical Solution of Burger's Huxley equation
%E*Uxx=Ut+ALPHA*UUx+BETA*(U^3+DELTA*U^2+GAMMA*U)
% , U(x,t)=GAMMA*[1+tanh(a1*(x-a2*t))]/2, t>=0
%--------------------------------------------------------------------------
clear all; % clear all variables in memory
LAM=1.6;
%LL=10;
%JJ=1000; % Solution at t=1
LL=8;
JJ=40;
H=1/LL;
K=LAM*H*H;
K=0.01;
%LLL=(LL/10)+9;
H= (2*pi)/100;
SGM=1.0;
EPS=1.0;
ALPHA=1.0;
BETA=0.0;
GAMMA=0.0;
DELTA=1.0;
a= (-1/4)*(1-SGM+(SGM)^2);
b= -SGM/(3*(1+SGM));
c= (-1/(2*SGM*(1+SGM)))*(1-SGM+(SGM)^2);
%a1=GAMMA*(-ALPHA+sqrt((ALPHA^2)+8*BETA))/8;
%a2=(ALPHA*GAMMA/2)+((2-GAMMA)*(ALPHA+sqrt((ALPHA^2)+8*BETA))/4);
%a1=(-ALPHA-sqrt((ALPHA^2)+8*BETA))/8;
%a2=(ALPHA/2)-((1-2*GAMMA)*(ALPHA-sqrt((ALPHA^2)+8*BETA))/4);
% LLL=(0.5*LL)+1;
%tt=(8.986818916*pi)/180;
%tt=8.986818916;
% %Spline in Tension
% tt=0.001;
% AA=(SGM*((tt/sin(tt))-cos(tt)))/(2*tt^2);
% BB1=(SGM*(-(tt*cos(tt)/sin(tt))+1))/(2*tt^2);
% BB2=(-(tt*cos(tt)/sin(tt))+1)/(2*tt^2);
% CC=((tt/sin(tt))-cos(tt))/(2*tt^2);
% %spline in compression
tt=(8.986818916*pi)/180;
AA=(SGM*((tt/sin(tt))-cos(tt)))/(2*tt^2);
BB1=(SGM*(-(tt*cos(tt)/sin(tt))+1))/(2*tt^2);
BB2=(-(tt*cos(tt)/sin(tt))+1)/(2*tt^2);
CC=((tt/sin(tt))-cos(tt))/(2*tt^2);
%Cubic Spline
%AA=SGM/3;
%BB1=SGM/6;
%BB2=1/6;
%CC=1/3;
disp('Numerical Solution of 1D Parabolic Equation');
disp(' ');
fprintf('LL=%4.2f, K=%6.4f, LAM=%6.2f \n',LL,K,LAM);
disp(' ');
% Dimensions
X=zeros(1,LL+1);
T=zeros(1,JJ+1);
TT=zeros(1,JJ+1);
U=zeros(LL+1,JJ+1);
EU=zeros(LL+1,JJ+1);
OLDU=zeros(LL+1,JJ+1);
H=zeros(1,LL);
if SGM==1
S=1;
for i=1:LL-1
S=S+SGM^i;
end;
% H(1)=2/S;
H(1)=1/S;
else
% H(1)=2*((1-SGM)/(1-(SGM^LL)));
H(1)=((1-SGM)/(1-(SGM^LL)));
end
% Formation of Variable step lengths
for L=1:LL-1
H(L+1)=SGM*H(L);
end
% Formation of Nodal Points
% X(1)=-1.0;
X(1)=0.0;
for L=2:LL
for J=1:JJ+1
X(L)=X(L-1)+H(L-1);
T(J)=(J-1)*K;
end
end
X(LL+1)=1.0;
for L=1:LL+1
for J=1:JJ+1
% EU(L,J)=GAMMA*(1+tanh(a1*(X(L)-a2*T(J))))/2;
% EU(L,J)=(1+tanh(a1*(X(L)-a2*T(J))))/2;
EU(L,J)= exp(-T(J))*sin(X(L));
end
end
% Initial conditions
for L=1:LL+1
%U(L,1)=EU(L,1);
U(L,1)=cos(pi*X(L));
end
% Boundary conditions
for J=2:JJ+1
%U(0,J)==EU(0,J);
U(1,J)=exp(-EPS*(pi)^2*T(J));
U(2,J)=-exp(-EPS*(pi)^2*T(J));
%U(1,J)=EU(1,J)==-exp(-EPS*(pi)^2*t);;
U(LL+1,J)=EU(LL+1,J);
end
% Solution at advanced time level
tic % CPU time starts
for J=1:JJ
% Initial & First Approx
for L=2:LL
OLDU(L,J+1)=EU(L,J+1);
U(L,J+1)=OLDU(L,J+1);
end
for KK=1:1200 %Iteration starts
% Solution at Internal grid points at each time level
for L=2:LL
M1 = U(L+1,J+1);
M2 = U(L+1,J);
M3 = U(L,J+1);
M4 = U(L,J);
M5 = U(L-1,J+1);
M6 = U(L-1,J);
N1 = (1/(2*K))*(M5+M3-M6-M4);
N2 = (M5+M3+M6+M4)*(M5+M3+M6+M4-2)*(M5+M3+M6+M4-2*GAMMA);
N3 = (M5+M3++M6+M4);
N4 = (1/(2*SGM*H(L)))*(M5+M6-M3-M4);
N5 = (1/K)*(M3-M4)+(BETA/8)*(M3+M4)*(M3+M4-2)*(M3+M4-2*GAMMA)+(ALPHA/2)*(M3+M4)*(1/(2*SGM*H(L)*(SGM+1)))*(M5+M6-(1-(SGM)^2)*(M3+M4)-((SGM)^2)*(M1+M2));
N6 = (1/(2*K))*(M5+M3-M6-M4)+(ALPHA/8)*(M5+M6+M3+M4)*(M5+M6+M3+M4)+(BETA/64)*(M5+M6+M3+M4)*(M5+M6+M3+M4-4)*(M5+M6+M3+M4-4*GAMMA);
N7 = (1/(2*K))*(M1+M3-M2-M4);
N8 = (M1+M3+M2+M4)*(M1+M3+M2+M4-2)*(M1+M3+M2+M4-2*GAMMA);
N9 = (M1+M3++M2+M4);
N10 = (1/(2*SGM*H(L)))*(M1+M2-M3-M4);
N11 = (1/K)*(M3-M4)+(BETA/8)*(M3+M4)*(M3+M4-2)*(M3+M4-2*GAMMA)+(ALPHA/2)*(M3+M4)*(1/(2*SGM*H(L)*(SGM+1)))*(M1+M2-(1-(SGM)^2)*(M3+M4)-((SGM)^2)*(M5+M6));
N12 = (1/(2*K))*(M1+M3-M2-M4)+(ALPHA/8)*(M1+M2+M3+M4)*(M1+M2+M3+M4)+(BETA/64)*(M1+M2+M3+M4)*(M1+M2+M3+M4-4)*(M1+M2+M3+M4-4*GAMMA);
N13 = (1/K)*(M3-M4);
N14 = (1/K)*(M1-M2)-(1+SGM)*N13+(SGM/K)*(M5-M6);
N15= (M3+M4)/(2);
N16 = M1+M2+(1+SGM)*(M3+M4)+SGM*(M5+M6);
N17 = M1+M2-(1-(SGM^2))*(M3+M4)-(SGM)^2*(M5+M6);
N18 = (1/(2*K))*(M1+M3-M2-M4);
N19 =N9^2;
N21 = M5+M3-M6-M4;
N22 = N3^2;
O1 = N1+(BETA/16)*(N2)+(ALPHA/4)*(N3)*(N4+(H(L)/4)*(2*BB2*N5-CC*N6));
O2 = N7+(BETA/16)*(N8)+(ALPHA/4)*(N9)*(N10+(H(L)/4)*(2*BB1*N11-AA*N12));
O3 = N13+c*N14+ALPHA*(N15+(a/SGM*(1+SGM))*N16)*(1/2*H(L)*SGM*(1+SGM))*N17 + b*H(L)*[N18+(ALPHA/8)*N19 + (BETA/16)*N9*(N9-4)*(N9-4*GAMMA)-(N21+(ALPHA/8)*(N3)^2)+(BETA/16)*(N3)*(N3-4)*(N3-4*GAMMA)]+BETA*[(N15+(a/(SGM*(SGM+1)))*N16)*(N15+(a/(2*SGM*(SGM+1)))*N16-1)*(N15+(a/(2*SGM*(SGM+1)))*N16-GAMMA)];
NAA = (EPS/2)*(M1+M2-2*M3-2*M4+M5+M6) - ((H.^2)/3)* (O1+O2+O3);

Risposte (2)

SOUMNATH PAUL
SOUMNATH PAUL il 7 Mag 2024
To differentiate (NAA) with respect to (U(L,J+1)) (denoted as (M3) in the given context), we need to apply the partial derivative to each term in the expression that contains (M3). When focusing on the direct appearances of (M3):
The term (\left(\frac{\text{EPS}}{2}\right) \cdot (M1 + M2 - 2 \cdot M3 - 2 \cdot M4 + M5 + M6)) directly includes (M3). The derivative of (-2 \cdot M3) with respect to (M3) is (-2).
For the terms within (O1), (O2), and (O3), we need to identify where (M3) appears and calculate the derivative accordingly. This part is complex due to the nonlinear nature of these expressions. However, since the question focuses on the differentiation with respect to (M3), we can simplify our task by noting that (M3) appears in many places within (N1), (N2), (N3), etc., contributing to (O1), (O2), and (O3).
Given the complexity and the nonlinear nature of (O1), (O2), and (O3), the partial derivative with respect to (M3) would involve taking the derivative of each of these expressions with respect to (M3), considering the chain rule and product rule where applicable.
Without the explicit forms of (O1), (O2), and (O3) simplified, the general approach to differentiate (NAA) with respect to (U(L,J+1)) involves:
Differentiating the linear part directly involving (M3), which gives (-2).
For each appearance of (M3) in (O1), (O2), and (O3), apply the derivative rules (chain rule, product rule) to find their contributions.
This approach outlines the methodology rather than providing a final numerical answer due to the complexity of the expressions involved. For a precise result, each term in (O1), (O2), and (O3) containing (M3) must be differentiated individually, and this process can be quite involved.
Hope it helps!
Regards,
Soumnath

John
John il 10 Mag 2024
To differentiate the non-linear advection-diffusion equation (NAA) with respect to U(L, J+1), we need to take the partial derivative of the equation with respect to U(L, J+1). The general form of the NAA equation is:
E * Uxx = Ut + ALPHA * U * Ux + BETA * (U^3 + DELTA * U^2 + GAMMA * U)
Here, E is the diffusion coefficient, ALPHA is the advection coefficient, and BETA, DELTA, and GAMMA are the non-linear coefficients.
To take the partial derivative with respect to U(L, J+1), we need to treat all other terms as constants. The derivative will be:
(E * Uxx)/U(L, J+1) = 0
(Ut)/U(L, J+1) = 1 (assuming an explicit time discretization)
(ALPHA * U * Ux)/U(L, J+1) = ALPHA * Ux
(BETA * (U^3 + DELTA * U^2 + GAMMA * U))/U(L, J+1) = BETA * (3 * U^2 + 2 * DELTA * U + GAMMA)
Combining these terms, we get:
(NAA)/U(L, J+1) = 1 + ALPHA * Ux + BETA * (3 * U^2 + 2 * DELTA * U + GAMMA)
To implement this in your code, you would need to calculate the spatial derivative Ux at the point (L, J+1) using a suitable finite difference approximation. Then, you can substitute the values of U(L, J+1), Ux(L, J+1), and the coefficients ALPHA, BETA, DELTA, and GAMMA into the above equation to obtain the desired derivative.
Note that this derivative would be used in the numerical solution of the NAA equation, typically in an iterative scheme or a Newton-like method to update the solution at each time step.

Categorie

Scopri di più su Partial Differential Equation Toolbox in Help Center e File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by