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
Version imprimable
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
Salut,
On aura du mal à te répondre. Montre au moins le fichier Mathematica.
(* 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
Quel est précisément le problème ?
je sais pas comment convertir une discrétisation en fortrant est ce que vous pouvez m'aider
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éou retardéCode:(F(x+dx)-F(x))/dx
d'ordre 1 ou plutôt centré d'ordre 2 (encore plus efficace)Code:(F(x+dx)-F(x))/dh
.Code:(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, après il y a beaucoup de schémas explicite et implicite .Code:(G(t+dt)-G(t))/dt
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....?
Citation:
Envoyé par IlhemADEL