Discussion :

[sayed]

j'ai un programme pour les equations d'euler j'arrive pas a l'exécuter!
comment declarer les variables d'entrées? je suis débutant en fortran ,
voici le code :

subroutine EULER2D(nx,ny,dx,dy,cfl,gamma,theta,tf,u)

Code :

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

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
c ******************************************************************
c * INPUT nx, ny: # of cells in x-, y-direction
c * dx, dy: step sizes in x-, y-direction
c * cfl: CFL # gamma: adiabatic constant of gas
c * tf: final time
c * theta=1: MM1 limiter; =2: MM2 limiter; >2: UNO limiter.
c * u: initial cell averages of conservative variables.
c * Supply entries of u((md+1):(nx+md),(md+1):(ny+md),4)
c * OUTPUT u: cell averages at final time "tf"
c * REMARK 1. Reset "nxd","nyd" to adjust array dimensions.
c * 2. Padded to each side of the computational domain are
c * "md" ghost cells, average values on which are
c * assigned by boundary conditions.
c ******************************************************************
c

parameter(md=3,mn=4,nxd=400+2*md,nyd=400+2*md)

Code :

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

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
c
      real u(nxd,nyd,mn), ux(nxd,nyd,mn), uy(nxd,nyd,mn)
      real f(nxd,nyd,mn), fx(nxd,nyd,mn)
      real g(nxd,nyd,mn), gy(nxd,nyd,mn)
      real v(nxd,nyd), du(nxd,2), df(nxd,2)
 
      xmin(a,b) = 0.5*(sign(1.,a)+sign(1.,b))*min(abs(a),abs(b))
      xmic(z,a,b) = xmin( z*xmin(a,b), 0.5*(a+b) )
      gm1 = gamma - 1.0
      tc = 0.0
      istop = 0
      do 1000 nt = 1, 10000
       do 999 io = 0, 1
c* Periodic boundary condition in both x- & y-direction
        do 101 m = 1, mn
         do 100 i = 1, md
          do 100 j = md + 1, ny + md
           u(i, j,m) = u(nx+i,j,m)
           u(nx + md+i,j,m) = u(md+i,j,m)
100       continue
        do 101 j = 1, md
         do 101 i = 1, nx + 2*md
          u(i,j, m) = u(i,ny+j,m)
          u(i,ny + md+j,m) = u(i,md+j,m)
101     continue
c* Compute f & g and maximum wave speeds "em_x, em_y".
c* See (2.1) & (4.3).
        em_x = 1.e-15
        em_y = 1.e-15
        do 200 j = 1, ny + 2*md
         do 200 i = 1, nx + 2*md
          den = u(i,j,1)
          vex = u(i,j,2) / den
          vey = u(i,j,3) / den
          eng = u(i,j,4)
          pre = gm1 * ( eng - .5*den*( vex*vex + vey*vey ) )
          cvel = sqrt( gamma * pre / den )
          em_x = max( em_x, abs(vex) + cvel )
          em_y = max( em_y, abs(vey) + cvel )
          f(i,j,1) = den * vex
          f(i,j,2) = den * vex**2 + pre
          f(i,j,3) = den * vex * vey
          f(i,j,4) = vex * ( pre + eng )
          g(i,j,1) = den * vey
          g(i,j,2) = den * vex * vey
          g(i,j,3) = den * vey**2 + pre
          g(i,j,4) = vey * ( pre + eng )
200     continue
c* Compute numerical derivatives "ux", "uy", "fx", "gy".
c*cSee (3.1) & (4.1)
        do 330 m = 1, mn
         do 310 j = 3, ny + 2*md - 2
          do 301 i = 1, nx + 2*md - 1
           du(i,1) = u(i+1,j,m) - u(i,j,m)
301        df(i,1) = f(i+1,j,m) - f(i,j,m)
          do 302 i = 1, nx + 2*md - 2
           du(i,2) = du(i+1,1) - du(i,1)
302        df(i,2) = df(i+1,1) - df(i,1)
           if( theta .lt. 2.5 ) then
            do 303 i = 3, nx + 2*md - 2
             ux(i,j,m) = xmic( theta, du(i-1,1), du(i,1) )
303          fx(i,j,m) = xmic( theta, df(i-1,1), df(i,1) )
           else
            do 304 i = 3, nx + 2*md - 2
             ux(i,j,m)=xmin(du(i-1,1)+.5*xmin(du(i-2,2),du(i-1,2)),
     &        du(i, 1)-.5*xmin(du(i-1,2),du(i, 2)))
             fx(i,j,m)=xmin(df(i-1,1)+.5*xmin(df(i-2,2),df(i-1,2)),
     &        df(i, 1)-.5*xmin(df(i-1,2),df(i, 2)))
304         continue
           endif
310      continue
         do 320 i = 3, nx + 2*md - 2
          do 311 j = 1, ny + 2*md - 1
           du(j,1) = u(i,j+1,m) - u(i,j,m)
311        df(j,1) = g(i,j+1,m) - g(i,j,m)
          do 312 j = 1, ny + 2*md - 2
           du(j,2) = du(j+1,1) - du(j,1)
312        df(j,2) = df(j+1,1) - df(j,1)
           if( theta .lt. 2.5 ) then
            do 313 j = 3, ny + 2*md - 2
             uy(i,j,m) = xmic( theta, du(j-1,1), du(j,1) )
313          gy(i,j,m) = xmic( theta, df(j-1,1), df(j,1) )
           else
            do 314 j = 3, ny + 2*md - 2
             uy(i,j,m)=xmin(du(j-1,1)+.5*xmin(du(j-2,2),du(j-1,2)),
     &        du(j, 1)-.5*xmin(du(j-1,2),du(j, 2)))
             gy(i,j,m)=xmin(df(j-1,1)+.5*xmin(df(j-2,2),df(j-1,2)),
     &        df(j, 1)-.5*xmin(df(j-1,2),df(j, 2)))
314         continue
           endif
320      continue
330     continue
c* Compute time step size according to the input CFL #.
        if(io.eq.0) then
         dt = cfl / max( em_x/dx, em_y/dy )
         if( ( tc + 2.*dt ) .ge. tf ) then
          dt = 0.5 * ( tf - tc )
          istop = 1
         endif
        endif
        dtcdx2 = 0.5 * dt / dx
        dtcdy2 = 0.5 * dt / dy
c* Compute the flux values of f & g at half time step.
c* See (2.15) & (2.16).
        do 400 j = 3, ny + 2*md - 2
         do 400 i = 3, nx + 2*md - 2
          den = u(i,j,1) - dtcdx2*fx(i,j,1) - dtcdy2*gy(i,j,1)
          xmt = u(i,j,2) - dtcdx2*fx(i,j,2) - dtcdy2*gy(i,j,2)
          ymt = u(i,j,3) - dtcdx2*fx(i,j,3) - dtcdy2*gy(i,j,3)
          eng = u(i,j,4) - dtcdx2*fx(i,j,4) - dtcdy2*gy(i,j,4)
          pre = gm1 * ( eng - .5 * ( xmt*xmt + ymt*ymt ) / den )
          f(i,j,1) = xmt
          f(i,j,2) = xmt * xmt / den + pre
          f(i,j,3) = xmt * ymt / den
          f(i,j,4) = xmt / den * ( pre + eng )
          g(i,j,1) = ymt
          g(i,j,2) = xmt * ymt / den
          g(i,j,3) = ymt * ymt / den + pre
          g(i,j,4) = ymt / den * ( pre + eng )
400     continue
c* Compute the values of "u" at the next time level. See (2.16).
        do 510 m = 1, mn
         do 501 j = md + 1 - io, ny + md - io
          do 501 i = md + 1 - io, nx + md - io
           v(i,j) = 0.25 * ( u(i,j, m) + u(i+1,j, m)
     &      + u(i,j+1,m) + u(i+1,j+1,m) )
     &      + 0.0625 * ( ux(i,j, m) - ux(i+1,j, m)
     &      + ux(i,j+1,m) - ux(i+1,j+1,m)
     &      + uy(i,j, m) + uy(i+1,j, m)
     &      - uy(i,j+1,m) - uy(i+1,j+1,m) )
     &      + dtcdx2 * ( f(i,j, m) - f(i+1,j, m)
     &      + f(i,j+1,m) - f(i+1,j+1,m) )
     &      + dtcdy2 * ( g(i,j, m) + g(i+1,j, m)
     &      - g(i,j+1,m) - g(i+1,j+1,m) )
501      continue
         do 502 j = md + 1, ny + md
          do 502 i = md + 1, nx + md
502        u(i,j,m) = v(i-io,j-io)
510     continue
        tc = tc + dt
999    continue
       if(istop.eq.1 ) goto 1001
1000  continue
1001  return
      end

[Dramac]

Hum à priori pas comme ça
une ligne de déclaration se fait comme ça

Code :

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

1
2
Type, dimension(le cas échéant) :: les variables
real, dimension(nxd,nyd,mn), ux(nxd,nyd,mn), uy(nxd,nyd,mn) :: u, ux, uy

Une "forme" de tableau par ligne.

En esperant que la réponse ne vienne pas trop tard

[Sylvain Bergeron]

La forme classique :

type variable1(dimensions), variable2(dimensions)

est toujours valide, même en F2003...