1 : C CONJ3W SOURCE CHAT 05/01/12 22:17:28 5004
2 : c----------------------------------------------------------------------
3 : c GENERAL DESCRIPTION:
4 : c This subroutine provides the jacobians which are the derivatives
5 : c of the numerical flux function defined at the cell interface
6 : c with respect to the conservative variables of the left and right
7 : c cells (relative to the cell interface).
8 : c
9 : c EQUATIONS: 3D Euler equations of gas dynamics
10 : c
11 : c
12 : c REFERENCE: JCP, 129, 364-382 (1996)
13 : c " A Sequel to AUSM: AUSM+ ".
14 : c----------------------------------------------------------------------
15 : c INPUT:
16 : c
17 : c alpha -- parameter of the AUSM+ scheme in the Pressure function;
18 : c ( -3/4 <= alpha <= 3/16 ) (imposed as a parameter)
19 : c
20 : c beta -- parameter of the AUSM+ scheme in the Mach function;
21 : c ( -1/16 <= beta <= 1/2 ) (imposed as a parameter)
22 : c
23 : c wvec_l -- vector of the primitive variables (rho,ux,uy,uz,p) at the
24 : c left cell;
25 : c
26 : c wvec_r -- vector of the primitive variables (rho,ux,uy,uz,p) at the
27 : c right cell;
28 : c
29 : c nvect -- normal vector to the interface (3 components in 3D);
30 : c
31 : c tvec1 -- tangential vector to the interface;
32 : c
33 : c tvec2 -- tangential vector to the interface;
34 : c
35 : c ga -- ratio of the specific heats (assumed constant)
36 : c----------------------------------------------------------------------
37 : c
38 : c OUTPUT:
39 : c
40 : c jtl -- jakobian matrix 5 by 5 - derivatives of the numerical
41 : c flux function with respect to the conservative variables
42 : c from the left cell;
43 : c
44 : c jtr -- jakobian matrix 5 by 5 - derivatives of the numerical
45 : c flux function with respect to the conservative variables
46 : c from the right cell.
47 : c----------------------------------------------------------------------
48 : SUBROUTINE CONJ3W(JTL,JTR,WVEC_L,WVEC_R,NVECT,TVECT1,TVECT2,
49 : & ga)
50 : IMPLICIT INTEGER(I-N)
51 : real*8 wvec_l(5),wvec_r(5)
52 : real*8 nvect(3),tvect1(3),tvect2(3)
53 : real*8 jl(5,5),jr(5,5)
54 : real*8 wl(5,5),wr(5,5)
55 : real*8 jtl(5,5),jtr(5,5)
56 : real*8 alpha,beta
57 : real*8 ga, gm1
58 : real*8 n_x,n_y,n_z,t1_x,t1_y,t1_z,t2_x,t2_y,t2_z
59 : real*8 un_l,un_r
60 : real*8 ml,mr,Mplus,Mmin
61 : real*8 am
62 : real*8 rold_l,uold_l,vold_l,wold_l,pold_l,eold_l
63 : real*8 rold_r,uold_r,vold_r,wold_r,pold_r,eold_r
64 : real*8 Pplus,Pmin
65 : real*8 temp_l,temp_r,brac_l,brac_r
66 : real*8 aleft, arigh
67 : real*8 damr_l,damr_r,damu_l,damu_r
68 : real*8 damv_l,damv_r,damp_l,damp_r
69 : real*8 damw_l,damw_r
70 : real*8 dmlr_l,dmlr_r,dmlu_l,dmlu_r
71 : real*8 dmlv_l,dmlv_r,dmlp_l,dmlp_r
72 : real*8 dmlw_l,dmlw_r
73 : real*8 dmrr_l,dmrr_r,dmru_l,dmru_r
74 : real*8 dmrv_l,dmrv_r,dmrp_l,dmrp_r
75 : real*8 dmrw_l,dmrw_r
76 : real*8 dPpr_l,dPpr_r,dPpu_l,dPpu_r
77 : real*8 dPpv_l,dPpv_r,dPpp_l,dPpp_r
78 : real*8 dPpw_l,dPpw_r
79 : real*8 dPmr_l,dPmr_r,dPmu_l,dPmu_r
80 : real*8 dPmv_l,dPmv_r,dPmp_l,dPmp_r
81 : real*8 dPmw_l,dPmw_r
82 : real*8 dpir_l,dpir_r,dpiu_l,dpiu_r
83 : real*8 dpiv_l,dpiv_r,dpip_l,dpip_r
84 : real*8 dpiw_l,dpiw_r
85 : real*8 c11,c12,c13,c21,c22,c23,c31,c32,c33,det
86 : real*8 zc11,zc12,zc13,zc21,zc22,zc23,zc31,zc32,zc33
87 : integer i,j,k
88 : parameter(alpha = 0.1875D0, beta = 0.125D0)
89 : c-----------------------
90 : gm1=ga-1.0d0
91 : c-----------------------
92 : n_x=nvect(1)
93 : n_y=nvect(2)
94 : n_z=nvect(3)
95 : c-----------------------
96 : t1_x=tvect1(1)
97 : t1_y=tvect1(2)
98 : t1_z=tvect1(3)
99 : c-----------------------
100 : t2_x=tvect2(1)
101 : t2_y=tvect2(2)
102 : t2_z=tvect2(3)
103 : c--------------------------------------
104 : rold_l=wvec_l(1)
105 : uold_l=wvec_l(2)
106 : vold_l=wvec_l(3)
107 : wold_l=wvec_l(4)
108 : pold_l=wvec_l(5)
109 : c-----------------------
110 : rold_r=wvec_r(1)
111 : uold_r=wvec_r(2)
112 : vold_r=wvec_r(3)
113 : wold_r=wvec_r(4)
114 : pold_r=wvec_r(5)
115 : c------------------------------------------------------------------
116 : c Computation of the specific total energy on the left and right.
117 : c------------------------------------------------------------------
118 : eold_l=(uold_l*uold_l+vold_l*vold_l+wold_l*wold_l)/2.0d0
119 : eold_l=eold_l+pold_l/(gm1*rold_l)
120 : eold_r=(uold_r*uold_r+vold_r*vold_r+wold_r*wold_r)/2.0d0
121 : eold_r=eold_r+pold_r/(gm1*rold_r)
122 : c-------------------------------------------------------------------
123 : c Computation of the speed of sound and its derivatives;
124 : c numerical speed of sound at the interface is taken as an average
125 : c of the speeds of sounds of the neighbouring cells
126 : c---------------------------------------------------------------------
127 : aleft=dsqrt(ga*pold_l/rold_l)
128 : arigh=dsqrt(ga*pold_r/rold_r)
129 : am=0.5d0*(aleft+arigh)
130 : c--------------------------------------------------------------------
131 : damr_r=-arigh/(4.0d0*rold_r)
132 : damu_r=0.0d0
133 : damv_r=0.0d0
134 : damw_r=0.0d0
135 : damp_r=ga/(4.0d0*arigh*rold_r)
136 : c-----------------------
137 : damr_l=-aleft/(4.0d0*rold_l)
138 : damu_l=0.0d0
139 : damv_l=0.0d0
140 : damw_l=0.0d0
141 : damp_l=ga/(4.0d0*aleft*rold_l)
142 : c-------------------------------------------------------------------
143 : c Computing numerical Mach number and its derivatives,
144 : c see p.370, under (A1).
145 : c-------------------------------------------------------------------
146 : un_l=uold_l*n_x+vold_l*n_y+wold_l*n_z
147 : un_r=uold_r*n_x+vold_r*n_y+wold_r*n_z
148 : c----------------------------------------
149 : ml=un_l/am
150 : mr=un_r/am
151 : c-------------------------------------------------------------------
151 : c Mplus and Mmin are calligraphic lettes M+ and M- from the paper,
152 : c see (19a) and (19b), p.367.
153 : c-------------------------------------------------------------------
154 : if(dabs(ml) .ge. 1.0d0) then
155 : Mplus=(ml+dabs(ml))/2.0d0
156 : else
157 : Mplus=(ml+1.0d0)*(ml+1.0d0)/4.0d0
158 : Mplus=Mplus+beta*(ml*ml-1.0d0)*(ml*ml-1.0d0)
159 : endif
160 : c-----------------
161 : if(dabs(mr) .ge. 1.0d0) then
162 : Mmin=(mr-dabs(mr))/2.0d0
163 : else
164 : Mmin=-(mr-1.0d0)*(mr-1.0d0)/4.0d0
165 : Mmin=Mmin-beta*(mr*mr-1.0d0)*(mr*mr-1.0d0)
166 : endif
167 : c-------------------------------------------------------------------
168 : c Derivatives of ml and mr with respect to both sets of primitive
169 : c variables: left and right.
170 : c-------------------------------------------------------------------
171 : temp_l=-un_l/(am*am)
172 : temp_r=-un_r/(am*am)
173 : c--------
174 : dmlr_l=temp_l*damr_l
175 : dmlr_r=temp_l*damr_r
176 : dmrr_l=temp_r*damr_l
177 : dmrr_r=temp_r*damr_r
178 : c--------
179 : dmlu_l=n_x/am+temp_l*damu_l
180 : dmlu_r=temp_l*damu_r
181 : dmru_l=temp_r*damu_l
182 : dmru_r=n_x/am+temp_r*damu_r
183 : c--------
184 : dmlv_l=n_y/am+temp_l*damv_l
185 : dmlv_r=temp_l*damv_r
186 : dmrv_l=temp_r*damv_l
187 : dmrv_r=n_y/am+temp_r*damv_r
188 : c--------
189 : dmlw_l=n_z/am+temp_l*damw_l
190 : dmlw_r=temp_l*damw_r
191 : dmrw_l=temp_r*damw_l
192 : dmrw_r=n_z/am+temp_r*damw_r
193 : c--------
194 : dmlp_l=temp_l*damp_l
195 : dmlp_r=temp_l*damp_r
196 : dmrp_l=temp_r*damp_l
197 : dmrp_r=temp_r*damp_r
198 : c---------------------------------------------------------------
198 : c Computing the calligraphic P+ and P- with their derivatives
199 : c see (21a) & (21b) on p.368.
200 : c---------------------------------------------------------------
201 : if(ml .ge. 1.0d0) then
202 : Pplus = 1.0d0
203 : else
204 : if((ml .gt. -1.0d0) .and. (ml .lt. 1.0d0)) then
205 : Pplus=(ml+1.0d0)*(ml+1.0d0)*(2.0d0-ml)/4.0d0
206 : Pplus=Pplus+alpha*ml*(ml*ml-1.0d0)*(ml*ml-1.0d0)
207 : else
208 : Pplus = 0.0d0
209 : endif
210 : endif
211 : c---------------------------------------------------------------
212 : if(mr .ge. 1.0d0) then
213 : Pmin = 0.0d0
214 : else
215 : if((mr .gt. -1.0d0) .and. (mr .lt. 1.0d0)) then
216 : Pmin=(mr-1.0d0)*(mr-1.0d0)*(2.0d0+mr)/4.0d0
217 : Pmin=Pmin-alpha*mr*(mr*mr-1.0d0)*(mr*mr-1.0d0)
218 : else
219 : Pmin = 1.0d0
220 : endif
221 : endif
222 : c---------------------------------------------------------------
223 : brac_l=(ml+1.0d0)*(2.0d0-ml)/2.0d0-(ml+1.0d0)*(ml+1.0d0)/4.0d0
224 : brac_l=brac_l+alpha*(ml*ml-1.0d0)*(ml*ml-1.0d0)
225 : brac_l=brac_l+4.0d0*alpha*ml*ml*(ml*ml-1.0d0)
226 : c--------------
227 : brac_r=(mr-1.0d0)*(2.0d0+mr)/2.0d0+(mr-1.0d0)*(mr-1.0d0)/4.0d0
228 : brac_r=brac_r-alpha*(mr*mr-1.0d0)*(mr*mr-1.0d0)
229 : brac_r=brac_r-4.0d0*alpha*mr*mr*(mr*mr-1.0d0)
230 : c---------------------------------------------------------------
231 : if((ml .gt. -1.0d0) .and. (ml .lt. 1.0d0)) then
232 : dPpr_l=brac_l*dmlr_l
233 : dPpr_r=brac_l*dmlr_r
234 : c------------
235 : dPpu_l=brac_l*dmlu_l
236 : dPpu_r=brac_l*dmlu_r
237 : c------------
238 : dPpv_l=brac_l*dmlv_l
239 : dPpv_r=brac_l*dmlv_r
240 : c------------
241 : dPpw_l=brac_l*dmlw_l
242 : dPpw_r=brac_l*dmlw_r
243 : c------------
244 : dPpp_l=brac_l*dmlp_l
245 : dPpp_r=brac_l*dmlp_r
246 : c------------
247 : else
248 : dPpr_l=0.0d0
249 : dPpr_r=0.0d0
250 : c-----------
251 : dPpu_l=0.0d0
252 : dPpu_r=0.0d0
253 : c-----------
254 : dPpv_l=0.0d0
255 : dPpv_r=0.0d0
256 : c-----------
257 : dPpw_l=0.0d0
258 : dPpw_r=0.0d0
259 : c-----------
260 : dPpp_l=0.0d0
261 : dPpp_r=0.0d0
262 : c-----------
263 : endif
264 : c---------------------------------------------------------------
265 : if((mr .gt. -1.0d0) .and. (mr .lt. 1.0d0)) then
266 : dPmr_l=brac_r*dmrr_l
267 : dPmr_r=brac_r*dmrr_r
268 : c------------
269 : dPmu_l=brac_r*dmru_l
270 : dPmu_r=brac_r*dmru_r
271 : c------------
272 : dPmv_l=brac_r*dmrv_l
273 : dPmv_r=brac_r*dmrv_r
274 : c------------
275 : dPmw_l=brac_r*dmrw_l
276 : dPmw_r=brac_r*dmrw_r
277 : c------------
278 : dPmp_l=brac_r*dmrp_l
279 : dPmp_r=brac_r*dmrp_r
280 : c------------
281 : else
282 : dPmr_l=0.0d0
283 : dPmr_r=0.0d0
284 : c-----------
285 : dPmu_l=0.0d0
286 : dPmu_r=0.0d0
287 : c-----------
288 : dPmv_l=0.0d0
289 : dPmv_r=0.0d0
290 : c-----------
291 : dPmw_l=0.0d0
292 : dPmw_r=0.0d0
293 : c-----------
294 : dPmp_l=0.0d0
295 : dPmp_r=0.0d0
296 : c-----------
297 : endif
298 : c---------------------------------------------------------------------
299 : c computing pmid -- p_{1/2} and its derivatives, see (20b), p.367.
300 : c---------------------------------------------------------------------
301 : dpir_l=dPpr_l*pold_l+dPmr_l*pold_r
302 : dpiu_l=dPpu_l*pold_l+dPmu_l*pold_r
303 : dpiv_l=dPpv_l*pold_l+dPmv_l*pold_r
304 : dpiw_l=dPpw_l*pold_l+dPmw_l*pold_r
305 : dpip_l=dPpp_l*pold_l+Pplus+dPmp_l*pold_r
306 : c----------------------------
307 : dpir_r=dPpr_r*pold_l+dPmr_r*pold_r
308 : dpiu_r=dPpu_r*pold_l+dPmu_r*pold_r
309 : dpiv_r=dPpv_r*pold_l+dPmv_r*pold_r
310 : dpiw_r=dPpw_r*pold_l+dPmw_r*pold_r
311 : dpip_r=dPpp_r*pold_l+Pmin+dPmp_r*pold_r
312 : c---------------------------------------------------------------------
313 : c computing JACOBIAN as a derivative of the numerical flux function
314 : c with respect to the primitive variables.
315 : c Notation: jl(i,j) --- is the derivative of the i-component of the
316 : c flux function with respect to the j-component of the
317 : c vector of primitive variables of the left state.
318 : c jr(i,j) --- is the derivative of the i-component of the
319 : c flux function with respect to the j-component of the
320 : c vector of primitive variables of the right state.
321 : c---------------------------------------------------------------------
322 : jl(1,1)=0.0D0
323 : jl(1,2)=0.0D0
324 : jl(1,3)=0.0D0
325 : jl(1,4)=0.0D0
326 : jl(1,5)=0.0D0
327 : c------------------------------------
328 : jr(1,1)=0.0D0
329 : jr(1,2)=0.0D0
330 : jr(1,3)=0.0D0
331 : jr(1,4)=0.0D0
332 : jr(1,5)=0.0D0
333 : c------------------------------------
334 : c------------------------------------
335 : jl(2,1)=dpir_l*n_x
336 : c------------------------
337 : jl(2,2)=dpiu_l*n_x
338 : c------------------------
339 : jl(2,3)=dpiv_l*n_x
340 : c------------------------
341 : jl(2,4)=dpiw_l*n_x
342 : c------------------------
343 : jl(2,5)=dpip_l*n_x
344 : c-------------------------------------------------
345 : c-------------------------------------------------
346 : jl(3,1)=dpir_l*n_y
347 : c------------------------
348 : jl(3,2)=dpiu_l*n_y
349 : c------------------------
350 : jl(3,3)=dpiv_l*n_y
351 : c------------------------
352 : jl(3,4)=dpiw_l*n_y
353 : c------------------------
354 : jl(3,5)=dpip_l*n_y
355 : c-------------------------------------------------
356 : c-------------------------------------------------
357 : jl(4,1)=dpir_l*n_z
358 : c------------------------
359 : jl(4,2)=dpiu_l*n_z
360 : c------------------------
361 : jl(4,3)=dpiv_l*n_z
362 : c------------------------
363 : jl(4,4)=dpiw_l*n_z
364 : c------------------------
365 : jl(4,5)=dpip_l*n_z
366 : c-------------------------------------------------------
367 : c derivatives with respect to the right
368 : c set of primitive variables
369 : c-------------------------------------------------------
370 : jr(2,1)=dpir_r*n_x
371 : c------------------------
372 : jr(2,2)=dpiu_r*n_x
373 : c------------------------
374 : jr(2,3)=dpiv_r*n_x
375 : c------------------------
376 : jr(2,4)=dpiw_r*n_x
377 : c------------------------
378 : jr(2,5)=dpip_r*n_x
379 : c-------------------------------------------------------
380 : jr(3,1)=dpir_r*n_y
381 : c------------------------
382 : jr(3,2)=dpiu_r*n_y
383 : c------------------------
384 : jr(3,3)=dpiv_r*n_y
385 : c------------------------
386 : jr(3,4)=dpiw_r*n_y
387 : c------------------------
388 : jr(3,5)=dpip_r*n_y
389 : c--------------------------------------------------------
390 : jr(4,1)=dpir_r*n_z
391 : c------------------------
392 : jr(4,2)=dpiu_r*n_z
393 : c------------------------
394 : jr(4,3)=dpiv_r*n_z
395 : c------------------------
396 : jr(4,4)=dpiw_r*n_z
397 : c------------------------
398 : jr(4,5)=dpip_r*n_z
399 : c-------------------------------------------------------
400 : c-------------------------------------------------------
401 : jl(5,1)=0.0D0
402 : jl(5,2)=0.0D0
403 : jl(5,3)=0.0D0
404 : jl(5,4)=0.0D0
405 : jl(5,5)=0.0D0
406 : c-------------------------------------------------
407 : c-------------------------------------------------
408 : jr(5,1)=0.0D0
409 : jr(5,2)=0.0D0
410 : jr(5,3)=0.0D0
411 : jr(5,4)=0.0D0
412 : jr(5,5)=0.0D0
413 : c-------------------------------------------------------------
414 : c matrix wl(i,j) represents the derivative of the i-component
415 : c of the vector of primitive variables of the left state with
416 : c respect to the j-component of the vector of the conservative
417 : c variables of the left state.
418 : c
419 : c Here: (rho, u, v, w, p) - vector of primitive variables;
420 : c (rho, rho u, rh o v, rho w, rho e) - conservative variables.
421 : c-------------------------------------------------------------
422 : wl(1,1)=1.0d0
423 : wl(1,2)=0.0d0
424 : wl(1,3)=0.0d0
425 : wl(1,4)=0.0d0
426 : wl(1,5)=0.0d0
427 : c------------------------------
428 : wl(2,1)=-uold_l/rold_l
429 : wl(2,2)=1.0d0/rold_l
430 : wl(2,3)=0.0d0
431 : wl(2,4)=0.0d0
432 : wl(2,5)=0.0d0
433 : c------------------------------
434 : wl(3,1)=-vold_l/rold_l
435 : wl(3,2)=0.0d0
436 : wl(3,3)=1.0d0/rold_l
437 : wl(3,4)=0.0d0
438 : wl(3,5)=0.0d0
439 : c------------------------------
440 : wl(4,1)=-wold_l/rold_l
441 : wl(4,2)=0.0d0
442 : wl(4,3)=0.0d0
443 : wl(4,4)=1.0d0/rold_l
444 : wl(4,5)=0.0d0
445 : c------------------------------
446 : wl(5,1)=gm1*(uold_l*uold_l+vold_l*vold_l+wold_l*wold_l)/2.0d0
447 : wl(5,2)=-uold_l*gm1
448 : wl(5,3)=-vold_l*gm1
449 : wl(5,4)=-wold_l*gm1
450 : wl(5,5)=gm1
451 : c------------------------------
452 : c------------------------------
453 : wr(1,1)=1.0d0
454 : wr(1,2)=0.0d0
455 : wr(1,3)=0.0d0
456 : wr(1,4)=0.0d0
457 : wr(1,5)=0.0d0
458 : c------------------------------
459 : wr(2,1)=-uold_r/rold_r
460 : wr(2,2)=1.0d0/rold_r
461 : wr(2,3)=0.0d0
462 : wr(2,4)=0.0d0
463 : wr(2,5)=0.0d0
464 : c------------------------------
465 : wr(3,1)=-vold_r/rold_r
466 : wr(3,2)=0.0d0
467 : wr(3,3)=1.0d0/rold_r
468 : wr(3,4)=0.0d0
469 : wr(3,5)=0.0d0
470 : c------------------------------
471 : wr(4,1)=-wold_r/rold_r
472 : wr(4,2)=0.0d0
473 : wr(4,3)=0.0d0
474 : wr(4,4)=1.0d0/rold_r
475 : wr(4,5)=0.0d0
476 : c------------------------------
477 : wr(5,1)=gm1*(uold_r*uold_r+vold_r*vold_r+wold_r*wold_r)/2.0d0
478 : wr(5,2)=-uold_r*gm1
479 : wr(5,3)=-vold_r*gm1
480 : wr(5,4)=-wold_r*gm1
481 : wr(5,5)=gm1
482 : c----------------------------------------------
483 : c----------------------------------------------
484 : do 1 i=1,5
485 : do 2 j=1,5
486 : jtl(i,j)=0.0d0
487 : jtr(i,j)=0.0d0
488 : do 3 k=1,5
489 : jtl(i,j)=jtl(i,j)+jl(i,k)*wl(k,j)
490 : jtr(i,j)=jtr(i,j)+jr(i,k)*wr(k,j)
491 : 3 continue
492 : 2 continue
493 : 1 continue
494 : c-----------------------------------------------
495 : c Taking in account the dependancy of variables
496 : c-----------------------------------------------
497 : c11=t1_y*t2_z - t1_z*t2_y
498 : c12=n_y*t2_z - n_z*t2_y
499 : c13=n_y*t1_z - n_z*t1_y
500 : c-------------------------------------
501 : c21=t1_x*t2_z - t1_z*t2_x
502 : c22=n_x*t2_z - n_z*t2_x
503 : c23=n_x*t1_z - n_z*t1_x
504 : c-------------------------------------
505 : c31=t1_x*t2_y - t1_y*t2_x
506 : c32=n_x*t2_y - n_y*t2_x
507 : c33=n_x*t1_y - n_y*t1_x
508 : det=n_x*c11 - n_y*c21 + n_z*c31
509 : c----------------------------------------------------------------------
510 : ZC11=-NVECT(1)*C11-TVECT1(1)*C12+TVECT2(1)*C13
511 : ZC12=-NVECT(2)*C11-TVECT1(2)*C12+TVECT2(2)*C13
512 : ZC13=-NVECT(3)*C11-TVECT1(3)*C12+TVECT2(3)*C13
513 : C---------------------------------
514 : ZC21=NVECT(1)*C21+TVECT1(1)*C22-TVECT2(1)*C23
515 : ZC22=NVECT(2)*C21+TVECT1(2)*C22-TVECT2(2)*C23
516 : ZC23=NVECT(3)*C21+TVECT1(3)*C22-TVECT2(3)*C23
517 : C---------------------------------
518 : ZC31=-NVECT(1)*C31-TVECT1(1)*C32+TVECT2(1)*C33
519 : ZC32=-NVECT(2)*C31-TVECT1(2)*C32+TVECT2(2)*C33
520 : ZC33=-NVECT(3)*C31-TVECT1(3)*C32+TVECT2(3)*C33
521 : c---------------------------------------------------------------------
522 : do 11 i=1,5
523 : jtl(i,1)=jtl(i,1)+jtr(i,1)
524 : jtl(i,2)=jtl(i,2)+jtr(i,2)*zc11/det+jtr(i,3)*zc21/det+
525 : & jtr(i,4)*zc31/det
526 : jtl(i,3)=jtl(i,3)+jtr(i,2)*zc12/det+jtr(i,3)*zc22/det+
527 : & jtr(i,4)*zc32/det
528 : jtl(i,4)=jtl(i,4)+jtr(i,2)*zc13/det+jtr(i,3)*zc23/det+
529 : & jtr(i,4)*zc33/det
530 : jtl(i,5)=jtl(i,5)+jtr(i,5)
531 : 11 continue
532 : c----------------------------------------------------------------------
533 : return
534 : end
535 :
536 :
537 :
538 :
539 :
540 :
541 :
542 :
543 :
544 :
545 :
546 :
547 :
548 :
549 :
550 :
551 :
552 :
553 :
© Cast3M 2003 - All rights reserved.
Disclaimer