Discussion :

[IlhemADEL]

bonjour j'ai un programme en mathématica et je dois le convertir en fortrant . du coup j'ai un problème du programmer une discrétisation sur le fortrant est ce vous avez une idée ?
merci

[__dardanos__]

Salut,
On aura du mal à te répondre. Montre au moins le fichier Mathematica.

[IlhemADEL]

(* DONNEES GENERALES *)

name = "ColE #9";
zmax = 250; (* Profondeur maximale de la simulation \
(mm)*)
tend = 3; (* temps total de la \
simulation (h)*)
zcol = 246.7; (* Profondeur de la colonne \
(mm)*)
ts = 1; (* durée de l'infiltration (h)*)
h = 25; (* pas d'espace (mm)*)
\[Tau] = 0.005; (* pas de temps (h)*)

S = 16733.06; (* surface colonne (mm^2)*)
puissance = 2;

(* DONNEES CONCERNANT LE MODELE KDW EAU *)

u0 = 0.001; (* intensité du flux \
d'eau initial dans la colonne (mm/h)*)
us = 19.5; (* intensité du flux d'eau imposé à \
la surface (mm/h)*)

a = 1.8; (*1.55*)
b = 10^4.5; (*10^3.8*)
\[Nu] = 40; (*40*)


p1 = b^(1/a); p2 = a b^(1/a); q = (a - 1)/a;

(* CALCUL DES VECTEURS TEMPS ET ESPACE *)

time = Table[j, {j, 0, tend, \[Tau]}];
z = Table[i, {i, 0, zmax, h}];
test1 = Abs[ts - time];
tf = First[Flatten[Position[test1, Min[test1]]]];
test2 = Abs[zcol - z];
zf = First[Flatten[Position[test2, Min[test2]]]];
ZM = Length[z];
TM = Length[time];


(* DECLARATION VECTEURS ET MATRICES *)
(* concernant le flux u[i,j] *) \
\
\
\
\
\
\
\
\
\
\
\
\
\
\


U = Array[u, {ZM + 1, TM + 1}, {0, 0}];
P0 = Array[op0, {ZM + 1}, {0}];

;
TP1 = Array[temp1, {ZM + 1}, {0}];

(* CONDITIONS INITIALES ET AUX LIMITES *)
(* Sur le flux u[i,j] *)
Do[u[i, 0] = u0, {i, 0, ZM + 1}];
Do[Which[0 <= j <= tf, u[0, j] = us, tf < j <= TM ,
u[0, j] = u0], {j, 0, TM}];

(* DEFINITIONS *)

speed[x_] := p (x)^q;


(* CALCUL DU FLUX D'EAU u[i,j] *)


Do[
(*Print["time= ",j];*)
Do[temp0[i] = u[i, j], {i, 0, ZM}];
Do[ If[(temp0[i] < us) && (j > tf), p = p2, p = p1];
c0[i] = speed[ temp0[i]], {i, 0, ZM}];

diag[1] = 1;(*diag*)
surdiag[1] = 0; (*surdiag*)

Do[
diag[i] = 1 + (2 \[Nu] c0[i - 1] \[Tau])/h^2;(*diag*)

sousdiag[i - 1] =
c0[i - 1] \[Tau] (-(1/(2 h)) - \[Nu] /h^2);(*sousdiag*)

surdiag[i] = c0[i - 1] \[Tau] (1/(2 h) - \[Nu] /h^2);(*surdiag*)

, {i, 2, ZM}];

diag[ZM + 1] = 1;(*diag*)
sousdiag[ZM] = 0;(*sousdiag*)

DIAG = Table[diag[i], {i, 1, ZM + 1}];
SOUSDIAG = Table[sousdiag[i], {i, 1, ZM}];
SURDIAG = Table[surdiag[i], {i, 1, ZM}];
Print[DIAG];
matrice =
SparseArray[{Band[{1, 2}] -> SOUSDIAG, Band[{1, 1}] -> DIAG,
Band[{2, 1}] -> SURDIAG}, ZM + 1];

TP1 = LinearSolve[matrice, TP0];


ClearAll[DIAG, diag, SOUSDIAG, sousdiag, SURDIAG, surdiag];

Do[
u[i - 1, j + 1] = TP1[[i]]
, {i, 2, ZM}
];


u[ZM, j + 1] = u[ZM - 1, j + 1];

, {j, 0, TM - 1}];

(*
<< PlotLegends`;
<< Graphics`Graphics`;*)

(*$TextStyle={FontFamily->"Times",FontSize->16};*)


test3 = zcol - z;
Clear[up];
Clear[uc];
Clear[ud];

(*<< LinearAlgebra`MatrixManipulation`;*)
iup = First[Flatten[Position[Negative[test3], True]]];
Print["iup", iup];
idown = iup - 1;
Print["idown", idown];
zup = z[[iup]];
Print["Ziup", zup]
zd = z[[idown]];
Print["Zidown", zd]
(*up=TakeRows[U,{iup,iup}];
ud=TakeRows[U,{idown,idown}];*)
up = Take[U, {iup, iup}];
ud = Take[U, {idown, idown}];
uc = Flatten[
ud*((zup - zcol)/(zup - zd)) + up*((zcol - zd)/(zup - zd))];




dat = Import["/Users/emichel/Desktop/Colonne E#9.txt", "Table"] ;

(*<< LinearAlgebra`MatrixManipulation`;
<<Statistics`DescriptiveStatistics`;
ttu=TakeColumns[dat,{1}];
qq=TakeColumns[dat,{2}];
lttu=Length[ttu];
hyd=AppendRows[ttu,qq];
*)

(*
ttu=dat[[All,1]];
qq=dat[[All,2]];
lttu=Length[ttu];
hyd=Append[ttu,qq];
*)
(* Calcul du RMSE *)
(*
<< LinearAlgebra`MatrixManipulation`;
Do[

ttu1= Table[ttu[[j]],{i,1,TM}];

test=ttu1-time;

iap=First[Flatten[Position[Negative[test],True]]];
iav=iap-1;tap=time[[iap]];tav=time[[iav]];
uap=uc[[iap]];uav=uc[[iav]];

mat=(tav uav
tap uap

);
inte=Interpolation[mat,InterpolationOrder->1];
ues[j]=inte[ttu[[j]]];

,{j,1,lttu}];(*ues[1]=u0;*)

ues1=Table[ues[i],{i,1,lttu}];

Y=(Abs[ues1-qq])^puissance;

(*<<Statistics`DescriptiveStatistics`;*)

Minimum=Mean[Flatten[Y]]*100/Mean[Flatten[qq]];
*)

fg1 = ListPlot[dat, Prolog -> AbsolutePointSize[4],
AxesOrigin -> {0, 0},
AxesLabel -> {"Time (h)",
"u (mm.\!\(\*SuperscriptBox[\(h\), \(-1\)]\))"},
PlotLabel -> name, DisplayFunction -> Identity];

gg1 = Table[{time[[j]], uc[[j]]}, {j, 1, TM}];

fg2 = ListPlot[gg1, PlotJoined -> True,
PlotStyle -> {RGBColor[1, 0, 0], Thickness[0.002]},
AxesOrigin -> {0, 0},
AxesLabel -> {"Time (h)",
"u (mm.\!\(\*SuperscriptBox[\(h\), \(-1\)]\))"},
DisplayFunction -> Identity];

Print["MA = ", Minimum];
Print["a = ", a];
Print["b = ", b];
Print["\[Nu] = ", \[Nu]];

Show[fg1, fg2, DisplayFunction -> $DisplayFunction]
Null

[__dardanos__]

Quel est précisément le problème ?

[IlhemADEL]

je sais pas comment convertir une discrétisation en fortrant est ce que vous pouvez m'aider

[mgoumine]

Salut , je fais les discrétisations sur Fortran90 . c'est à dire -je sais pas si j'arriverai à vous aider- lorsque vous avez un opérateur différentiel tel que Grad , Rot, ou tout simplement une dérivée spatiale d'une fonction F(x) ...et ben celle ci par exemple se traduit en discrétisation par un schéma différentiel avancé

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(F(x+dx)-F(x))/dx

ou retardé

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(F(x+dx)-F(x))/dh

d'ordre 1 ou plutôt centré d'ordre 2 (encore plus efficace)

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(F(x+dx)-F(x-dx)/2dx

.

D'une autre part si vous avez aussi une dérivée temporaire d'une fonction quelconque encore qu'on appelle G(t) qui s'introduit dans l'équation à discrétiser alors elle devient

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(G(t+dt)-G(t))/dt

, après il y a beaucoup de schémas explicite et implicite .

En fin de compte vous appliquez ces schéma à l'équation de départ pour les différents temps (itérations) pour obtenir un système souvent matriciel à résoudre .

N.B: souvent le plus délicat c'est les formes que prennent cette équation de départ au instant de départ et aux bords .
N.B: vous allez discrétiser en utilisant quelle méthode SVP? Volume finis Éléments Finis....?

je sais pas comment convertir une discrétisation en fortrant est ce que vous pouvez m'aider