5 t7 =
rx(
i + 1,
j,0,0) ** 2
6 t8 =
rx(i + 1,j,0,1) ** 2
8 t13 =
rx(i,j,0,0) *
rx(i,j,1,1) -
rx(i,j,1,0) *
rx(i,j,0,1)
10 t15 =
rx(i,j,0,0) ** 2
11 t16 =
rx(i,j,0,1) ** 2
12 t18 = t14 * (t15 + t16)
14 t20 = t2 * t19 / 0.2D1
16 t23 = (
u(i + 1,j,
n) -
u(i,j,
n)) * t22
18 t25 =
ut(i + 1,j,
n) -
ut(i,j,
n)
19 t34 = 0.1D1 / (
rx(i + 2,j,0,0) *
rx(i + 2,j,1,1) -
rx(i + 2,j,1,0)
21 t35 =
rx(i + 2,j,0,0) ** 2
22 t36 =
rx(i + 2,j,0,1) ** 2
23 t42 = (
u(i + 2,j,
n) -
u(i + 1,j,
n)) * t22
26 t61 =
rx(i + 1,j,0,0) *
rx(i + 1,j,1,0) +
rx(i + 1,j,0,1) *
rx(i +
28 t63 = (
u(i + 1,j + 1,
n) -
u(i + 1,j,
n)) * t51
29 t65 = (
u(i + 1,j,
n) -
u(i + 1,j - 1,
n)) * t51
30 t68 = t58 * t61 * (t63 + t65) / 0.2D1
32 t75 =
rx(i,j,0,0) *
rx(i,j,1,0) +
rx(i,j,0,1) *
rx(i,j,1,1)
33 t77 = (
u(i,j + 1,
n) -
u(i,j,
n)) * t51
34 t79 = (
u(i,j,
n) -
u(i,j - 1,
n)) * t51
35 t82 = t72 * t75 * (t77 + t79) / 0.2D1
36 t85 = (t68 - t82) * t22 / 0.2D1
37 t89 = 0.1D1 / (
rx(i + 1,j + 1,0,0) *
rx(i + 1,j + 1,1,1) -
rx(i +
38 #1,j + 1,1,0) * rx(i + 1,j + 1,0,1))
40 t93 =
rx(i + 1,j + 1,0,0) *
rx(i + 1,j + 1,1,0) +
rx(i + 1,j + 1,0
41 #,1) *
rx(i + 1,j + 1,1,1)
42 t97 = (
u(i + 1,j + 1,
n) -
u(i,j + 1,
n)) * t22
43 t103 = t58 * t61 * (t42 + t23) / 0.2D1
44 t110 = 0.1D1 / (
rx(i + 1,j - 1,0,0) *
rx(i + 1,j - 1,1,1) -
rx(i +
45 # 1,j - 1,1,0) * rx(i + 1,j - 1,0,1))
47 t114 =
rx(i + 1,j - 1,0,0) *
rx(i + 1,j - 1,1,0) +
rx(i + 1,j - 1,
48 #0,1) *
rx(i + 1,j - 1,1,1)
49 t118 = (
u(i + 1,j - 1,
n) -
u(i,j - 1,
n)) * t22
50 t125 =
rx(i + 1,j + 1,1,0) ** 2
51 t126 =
rx(i + 1,j + 1,1,1) ** 2
52 t129 =
rx(i + 1,j,1,0) ** 2
53 t130 =
rx(i + 1,j,1,1) ** 2
54 t132 = t6 * (t129 + t130)
55 t136 =
rx(i + 1,j - 1,1,0) ** 2
56 t137 =
rx(i + 1,j - 1,1,1) ** 2
58 t157 =
rx(i - 1,j,0,0) *
rx(i - 1,j,1,1) -
rx(i - 1,j,1,0) *
rx(i
61 t159 =
rx(i - 1,j,0,0) ** 2
62 t160 =
rx(i - 1,j,0,1) ** 2
63 t162 = t158 * (t159 + t160)
65 t164 = t2 * t163 / 0.2D1
66 t166 = (
u(i,j,
n) -
u(i - 1,j,
n)) * t22
69 t173 =
rx(i - 1,j,0,0) *
rx(i - 1,j,1,0) +
rx(i - 1,j,0,1) *
rx(i
71 t175 = (
u(i - 1,j + 1,
n) -
u(i - 1,j,
n)) * t51
72 t177 = (
u(i - 1,j,
n) -
u(i - 1,j - 1,
n)) * t51
73 t180 = t170 * t173 * (t175 + t177) / 0.2D1
74 t183 = (t82 - t180) * t22 / 0.2D1
75 t186 =
rx(i,j + 1,0,0) *
rx(i,j + 1,1,1) -
rx(i,j + 1,1,0) *
rx(i,
79 t191 =
rx(i,j + 1,0,0) *
rx(i,j + 1,1,0) +
rx(i,j + 1,0,1) *
rx(i,
81 t193 = (
u(i,j + 1,
n) -
u(i - 1,j + 1,
n)) * t22
82 t196 = t188 * t191 * (t97 + t193) / 0.2D1
83 t199 = t72 * t75 * (t23 + t166) / 0.2D1
84 t202 = (t196 - t199) * t51 / 0.2D1
85 t205 =
rx(i,j - 1,0,0) *
rx(i,j - 1,1,1) -
rx(i,j - 1,1,0) *
rx(i,
89 t210 =
rx(i,j - 1,0,0) *
rx(i,j - 1,1,0) +
rx(i,j - 1,0,1) *
rx(i,
91 t212 = (
u(i,j - 1,
n) -
u(i - 1,j - 1,
n)) * t22
92 t215 = t207 * t210 * (t118 + t212) / 0.2D1
93 t218 = (t199 - t215) * t51 / 0.2D1
94 t219 =
rx(i,j + 1,1,0) ** 2
95 t220 =
rx(i,j + 1,1,1) ** 2
96 t222 = t187 * (t219 + t220)
97 t223 =
rx(i,j,1,0) ** 2
98 t224 =
rx(i,j,1,1) ** 2
99 t226 = t14 * (t223 + t224)
101 t228 = t2 * t227 / 0.2D1
103 t230 =
rx(i,j - 1,1,0) ** 2
104 t231 =
rx(i,j - 1,1,1) ** 2
105 t233 = t206 * (t230 + t231)
107 t235 = t2 * t234 / 0.2D1
109 t242 = dt * ((t24 - t167) * t22 + t85 + t183 + t202 + t218 + (t229
110 # - t236) * t51) * t13 / 0.2D1
111 t243 =
ut(i,j,
n) -
ut(i - 1,j,
n)
113 t247 =
dx * (t151 + t244) / 0.4D1
114 t254 = t24 + t20 * dt * t25 * t22 / 0.2D1 +
cc * t19 * (
ut(i + 1,j
115 #,n) + dt * ((t2 * (t34 * (t35 + t36) + t10) * t42 / 0.2D1 - t24) *
116 # t22 + (t2 * t34 * (rx(i + 2,j,0,0) * rx(i + 2,j,1,0) + rx(i + 2,j
117 #,0,1) * rx(i + 2,j,1,1)) * ((u(i + 2,j + 1,n) - u(i + 2,j,n)) * t5
118 #1 + (u(i + 2,j,n) - u(i + 2,j - 1,n)) * t51) / 0.2D1 - t68) * t22
119 #/ 0.2D1 + t85 + (t90 * t93 * ((u(i + 2,j + 1,n) - u(i + 1,j + 1,n)
120 #) * t22 + t97) / 0.2D1 - t103) * t51 / 0.2D1 + (t103 - t111 * t114
121 # * ((u(i + 2,j - 1,n) - u(i + 1,j - 1,n)) * t22 + t118) / 0.2D1) *
122 # t51 / 0.2D1 + (t2 * (t89 * (t125 + t126) + t132) * t63 / 0.2D1 -
123 #t2 * (t132 + t110 * (t136 + t137)) * t65 / 0.2D1) * t51) * t5 / 0.
124 #2D1 - dx * ((ut(i + 2,j,n) - ut(i + 1,j,n)) * t22 + t151) / 0.4D1
125 #- ut(i,j,n) - t242 - t247) * sqrt(0.2D1) * (t7 + t8 + t15 + t16) *
126 #* (-0.1D1 / 0.2D1) / 0.4D1
129 t260 = t2 * (t6 * t61 + t258) / 0.2D1
130 t267 =
ut(i,j + 1,
n) -
ut(i,j,
n)
132 t269 =
ut(i,j,
n) -
ut(i,j - 1,
n)
134 t275 = t260 * (t63 + t65 + t77 + t79) / 0.4D1 + t260 * dt * ((
ut(i
135 # + 1,j + 1,n) - ut(i + 1,j,n)) * t51 + (ut(i + 1,j,n) - ut(i + 1,j
136 # - 1,n)) * t51 + t268 + t270) / 0.8D1
137 t285 = 0.1D1 / (
rx(i - 2,j,0,0) *
rx(i - 2,j,1,1) -
rx(i - 2,j,1,0
138 #) *
rx(i - 2,j,0,1))
139 t286 =
rx(i - 2,j,0,0) ** 2
140 t287 =
rx(i - 2,j,0,1) ** 2
141 t293 = (
u(i - 1,j,
n) -
u(i - 2,j,
n)) * t22
142 t314 = 0.1D1 / (
rx(i - 1,j + 1,0,0) *
rx(i - 1,j + 1,1,1) -
rx(i -
143 # 1,j + 1,1,0) *
rx(i - 1,j + 1,0,1))
145 t318 =
rx(i - 1,j + 1,0,0) *
rx(i - 1,j + 1,1,0) +
rx(i - 1,j + 1,
146 #0,1) *
rx(i - 1,j + 1,1,1)
147 t326 = t170 * t173 * (t166 + t293) / 0.2D1
148 t333 = 0.1D1 / (
rx(i - 1,j - 1,0,0) *
rx(i - 1,j - 1,1,1) -
rx(i -
149 # 1,j - 1,1,0) * rx(i - 1,j - 1,0,1))
151 t337 =
rx(i - 1,j - 1,0,0) *
rx(i - 1,j - 1,1,0) +
rx(i - 1,j - 1,
152 #0,1) *
rx(i - 1,j - 1,1,1)
153 t346 =
rx(i - 1,j + 1,1,0) ** 2
154 t347 =
rx(i - 1,j + 1,1,1) ** 2
155 t350 =
rx(i - 1,j,1,0) ** 2
156 t351 =
rx(i - 1,j,1,1) ** 2
157 t353 = t158 * (t350 + t351)
158 t357 =
rx(i - 1,j - 1,1,0) ** 2
159 t358 =
rx(i - 1,j - 1,1,1) ** 2
160 t381 = t167 + t164 * dt * t243 * t22 / 0.2D1 +
cc * t163 * (
ut(i,j
161 #,n) + t242 - t247 - ut(i - 1,j,n) - dt * ((t167 - t2 * (t162 + t28
162 #5 * (t286 + t287)) * t293 / 0.2D1) * t22 + t183 + (t180 - t2 * t28
163 #5 * (rx(i - 2,j,0,0) * rx(i - 2,j,1,0) + rx(i - 2,j,0,1) * rx(i -
164 #2,j,1,1)) * ((u(i - 2,j + 1,n) - u(i - 2,j,n)) * t51 + (u(i - 2,j,
165 #n) - u(i - 2,j - 1,n)) * t51) / 0.2D1) * t22 / 0.2D1 + (t315 * t31
166 #8 * (t193 + (u(i - 1,j + 1,n) - u(i - 2,j + 1,n)) * t22) / 0.2D1 -
167 # t326) * t51 / 0.2D1 + (t326 - t334 * t337 * (t212 + (u(i - 1,j -
168 #1,n) - u(i - 2,j - 1,n)) * t22) / 0.2D1) * t51 / 0.2D1 + (t2 * (t3
169 #14 * (t346 + t347) + t353) * t175 / 0.2D1 - t2 * (t353 + t333 * (t
170 #357 + t358)) * t177 / 0.2D1) * t51) * t157 / 0.2D1 - dx * (t244 +
171 #(ut(i - 1,j,n) - ut(i - 2,j,n)) * t22) / 0.4D1) * sqrt(0.2D1) * (t
172 #15 + t16 + t159 + t160) ** (-0.1D1 / 0.2D1) / 0.4D1
173 t385 = t2 * (t258 + t158 * t173) / 0.2D1
174 t396 = t385 * (t77 + t79 + t175 + t177) / 0.4D1 + t385 * dt * (t26
175 #8 + t270 + (ut(i - 1,j + 1,n) - ut(i - 1,j,n)) * t51 + (ut(i - 1,j
176 #,n) - ut(i - 1,j - 1,n)) * t51) / 0.8D1
177 t403 = t2 * (t187 * t191 + t258) / 0.2D1
178 t414 = t403 * (t97 + t193 + t23 + t166) / 0.4D1 + t403 * dt * ((
ut
179 #(i + 1,j + 1,n) - ut(i,j + 1,n)) * t22 + (ut(i,j + 1,n) - ut(i - 1
180 #,j + 1,n)) * t22 + t151 + t244) / 0.8D1
181 t421 =
rx(i + 1,j + 1,0,0) ** 2
182 t422 =
rx(i + 1,j + 1,0,1) ** 2
183 t425 =
rx(i,j + 1,0,0) ** 2
184 t426 =
rx(i,j + 1,0,1) ** 2
185 t428 = t187 * (t425 + t426)
186 t432 =
rx(i - 1,j + 1,0,0) ** 2
187 t433 =
rx(i - 1,j + 1,0,1) ** 2
188 t447 = (
u(i,j + 2,
n) -
u(i,j + 1,
n)) * t51
189 t450 = t188 * t191 * (t447 + t77) / 0.2D1
190 t465 = 0.1D1 / (
rx(i,j + 2,0,0) *
rx(i,j + 2,1,1) -
rx(i,j + 2,1,0
191 #) *
rx(i,j + 2,0,1))
192 t480 =
rx(i,j + 2,1,0) ** 2
193 t481 =
rx(i,j + 2,1,1) ** 2
194 t500 = dy * (t268 + t270) / 0.4D1
195 t507 = t229 + t228 * dt * t267 * t51 / 0.2D1 +
cc * t227 * (
ut(i,j
196 # + 1,
n) + dt * ((t2 * (t89 * (t421 + t422) + t428) * t97 / 0.2D1 -
197 # t2 * (t428 + t314 * (t432 + t433)) * t193 / 0.2D1) * t22 + (t90 *
198 # t93 * ((u(i + 1,j + 2,n) - u(i + 1,j + 1,n)) * t51 + t63) / 0.2D1
199 # - t450) * t22 / 0.2D1 + (t450 - t315 * t318 * ((u(i - 1,j + 2,n)
200 #- u(i - 1,j + 1,n)) * t51 + t175) / 0.2D1) * t22 / 0.2D1 + (t2 * t
201 #465 * (rx(i,j + 2,0,0) * rx(i,j + 2,1,0) + rx(i,j + 2,0,1) * rx(i,
202 #j + 2,1,1)) * ((u(i + 1,j + 2,n) - u(i,j + 2,n)) * t22 + (u(i,j +
203 #2,n) - u(i - 1,j + 2,n)) * t22) / 0.2D1 - t196) * t51 / 0.2D1 + t2
204 #02 + (t2 * (t465 * (t480 + t481) + t222) * t447 / 0.2D1 - t229) *
205 #t51) * t186 / 0.2D1 - dy * ((ut(i,j + 2,n) - ut(i,j + 1,n)) * t51
206 #+ t268) / 0.4D1 - ut(i,j,n) - t242 - t500) * sqrt(0.2D1) * (t219 +
207 # t220 + t223 + t224) ** (-0.1D1 / 0.2D1) / 0.4D1
208 t511 = t2 * (t258 + t206 * t210) / 0.2D1
209 t522 = t511 * (t23 + t166 + t118 + t212) / 0.4D1 + t511 * dt * (t1
210 #51 + t244 + (ut(i + 1,j - 1,n) - ut(i,j - 1,n)) * t22 + (ut(i,j -
211 #1,n) - ut(i - 1,j - 1,n)) * t22) / 0.8D1
212 t529 =
rx(i + 1,j - 1,0,0) ** 2
213 t530 =
rx(i + 1,j - 1,0,1) ** 2
214 t533 =
rx(i,j - 1,0,0) ** 2
215 t534 =
rx(i,j - 1,0,1) ** 2
216 t536 = t206 * (t533 + t534)
217 t540 =
rx(i - 1,j - 1,0,0) ** 2
218 t541 =
rx(i - 1,j - 1,0,1) ** 2
219 t555 = (
u(i,j - 1,
n) -
u(i,j - 2,
n)) * t51
220 t558 = t207 * t210 * (t79 + t555) / 0.2D1
221 t573 = 0.1D1 / (
rx(i,j - 2,0,0) *
rx(i,j - 2,1,1) -
rx(i,j - 2,1,0
222 #) *
rx(i,j - 2,0,1))
223 t588 =
rx(i,j - 2,1,0) ** 2
224 t589 =
rx(i,j - 2,1,1) ** 2
225 t612 = t236 + t235 * dt * t269 * t51 / 0.2D1 +
cc * t234 * (
ut(i,j
226 #,
n) + t242 - t500 -
ut(i,j - 1,
n) - dt * ((t2 * (t110 * (t529 + t5
227 #30) + t536) * t118 / 0.2D1 - t2 * (t536 + t333 * (t540 + t541)) *
228 #t212 / 0.2D1) * t22 + (t111 * t114 * (t65 + (u(i + 1,j - 1,n) - u(
229 #i + 1,j - 2,n)) * t51) / 0.2D1 - t558) * t22 / 0.2D1 + (t558 - t33
230 #4 * t337 * (t177 + (u(i - 1,j - 1,n) - u(i - 1,j - 2,n)) * t51) /
231 #0.2D1) * t22 / 0.2D1 + t218 + (t215 - t2 * t573 * (rx(i,j - 2,0,0)
232 # * rx(i,j - 2,1,0) + rx(i,j - 2,0,1) * rx(i,j - 2,1,1)) * ((u(i +
233 #1,j - 2,n) - u(i,j - 2,n)) * t22 + (u(i,j - 2,n) - u(i - 1,j - 2,n
234 #)) * t22) / 0.2D1) * t51 / 0.2D1 + (t236 - t2 * (t233 + t573 * (t5
235 #88 + t589)) * t555 / 0.2D1) * t51) * t205 / 0.2D1 - dy * (t270 + (
236 #ut(i,j - 1,n) - ut(i,j - 2,n)) * t51) / 0.4D1) * sqrt(0.2D1) * (t2
237 #23 + t224 + t230 + t231) ** (-0.1D1 / 0.2D1) / 0.4D1
238 cg(1,1) =
u(i,j,
n) + dt *
ut(i,j,
n) + (t254 * t255 + t275 * t255 -
239 # t381 * t255 - t396 * t255) * t13 * t22 / 0.2D1 + (t414 * t255 + t
240 #507 * t255 - t522 * t255 - t612 * t255) * t13 * t51 / 0.2D1
241 cg(2,1) =
ut(i,j,
n) + (t254 * dt + t275 * dt - t381 * dt - t396 *
242 #dt) * t13 * t22 + (t414 * dt + t507 * dt - t522 * dt - t612 * dt)