SFEMaNS  version 4.1 (work in progress)
Reference documentation for SFEMaNS
 All Classes Files Functions Variables Groups Pages
post_processing_debug.f90
Go to the documentation of this file.
2 
3 PUBLIC :: post_proc_test
4 
5 PRIVATE
6 CONTAINS
7  !---------------------------------------------------------------------------
8  SUBROUTINE post_proc_test(vv_mesh, pp_mesh, temp_mesh, H_mesh, phi_mesh, list_mode, &
9  un, pn, hn, bn, phin, tempn, level_setn, mu_h_field, &
10  time, m_max_c, comm_one_d, comm_one_d_ns)
11  USE boundary
12  USE def_type_mesh
13  USE input_data
14  USE my_util
15  USE tn_axi
17  USE sft_parallele
18  IMPLICIT NONE
19  TYPE(mesh_type), POINTER :: pp_mesh, vv_mesh
20  TYPE(mesh_type), POINTER :: temp_mesh
21  TYPE(mesh_type), POINTER :: h_mesh, phi_mesh
22  INTEGER, POINTER, DIMENSION(:) :: list_mode
23  REAL(KIND=8), POINTER, DIMENSION(:,:,:) :: un, pn, hn, bn, phin, tempn
24  REAL(KIND=8), POINTER, DIMENSION(:,:,:,:) :: level_setn
25  REAL(KIND=8), POINTER, DIMENSION(:) :: mu_h_field
26  REAL(KIND=8) :: time
27  INTEGER :: m_max_c
28  REAL(KIND=8), DIMENSION(SIZE(un,1), SIZE(un,2), SIZE(un,3)) :: un_m1, un_ex, un_error
29  REAL(KIND=8), DIMENSION(SIZE(pn,1), SIZE(pn,2), SIZE(pn,3)) :: pn_m1, pn_ex, pn_error
30  REAL(KIND=8), DIMENSION(SIZE(Hn,1), SIZE(Hn,2), SIZE(Hn,3)) :: hn1, hn_ex, hn_error
31  REAL(KIND=8), DIMENSION(SIZE(phin,1),SIZE(phin,2),SIZE(phin,3)) :: phin1
32  REAL(KIND=8), DIMENSION(SIZE(tempn,1), SIZE(tempn,2), SIZE(tempn,3)) :: tempn_m1, tempn_ex, tempn_error
33  REAL(KIND=8), DIMENSION(SIZE(level_setn,1),SIZE(level_setn,2),SIZE(level_setn,3),SIZE(level_setn,4)):: level_setn_m1
34 
35  INTEGER :: i, k, code, rank, int_nb, n, type
36  REAL(KIND=8), DIMENSION(4) :: norm_err_loc, norm_err
37  REAL(KIND=8) :: err, norm
38  REAL(KIND=8) :: moyenne
39 !!$ REAL(KIND=8) :: error_on_tests=1.d-6
40 #include "petsc/finclude/petsc.h"
41  mpi_comm, DIMENSION(:), POINTER :: comm_one_d, comm_one_d_ns
42 
43  CALL mpi_comm_rank(mpi_comm_world,rank,code)
44 
45  SELECT CASE(inputs%numero_du_test_debug)
46  CASE(1,2,8,9,19,20,24,25)
47  DO i = 1, m_max_c
48  DO k= 1, 6
49  un_m1(:,k,i) = un(:,k,i) - vv_exact(k,vv_mesh%rr,list_mode(i),time)
50  END DO
51  DO k= 1, 2
52  pn_m1(:,k,i) = pn(:,k,i) - pp_exact(k,pp_mesh%rr,list_mode(i),time)
53  END DO
54  IF (list_mode(i) == 0) THEN
55  CALL moy(comm_one_d(1),pp_mesh, pn_m1(:,1,i),moyenne)
56  pn_m1(:,1,i) = pn_m1(:,1,i) - moyenne
57  ENDIF
58  END DO
59 
60  !norm_err(1) = norm_SF(comm_one_d_NS, 'L2', vv_mesh, list_mode, un_m1)
61  norm_err(1) = sqrt(dot_product_sf(comm_one_d_ns,vv_mesh, list_mode, un_m1, un_m1))
62  norm_err(2) = norm_sf(comm_one_d_ns, 'sH1', vv_mesh, list_mode, un_m1)
63  norm_err(3) = norm_sf(comm_one_d_ns, 'div', vv_mesh, list_mode, un)
64  norm_err(4) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, pn_m1)
65  IF (rank==0) THEN
66  WRITE(10,*) 'Velocity field #####################'
67  WRITE(10,*) 'L2 error on velocity = ', norm_err(1)
68  WRITE(10,*) 'H1 error on velocity = ', norm_err(2)
69  WRITE(10,*) 'L2 norm of divergence = ', norm_err(3)
70  WRITE(10,*) 'Pressure field #####################'
71  WRITE(10,*) 'L2 error on pressure = ', norm_err(4)
72  END IF
73 
74  IF (inputs%if_temperature) THEN
75  DO i = 1, m_max_c
76  DO k= 1, 2
77  tempn_m1(:,k,i) = tempn(:,k,i) - temperature_exact(k,temp_mesh%rr, list_mode(i), time)
78  END DO
79  END DO
80  norm_err(1) = norm_sf(comm_one_d_ns, 'L2', temp_mesh, list_mode, tempn_m1)
81  norm_err(2) = norm_sf(comm_one_d_ns, 'sH1', temp_mesh, list_mode, tempn_m1)
82  IF (rank==0) THEN
83  WRITE(10,*) 'Temperature field#####################'
84  WRITE(10,*) 'L2 error on temperature = ', norm_err(1)
85  WRITE(10,*) 'H1 error on temperature = ', norm_err(2)
86  WRITE(10,*) 'Temperature field#####################'
87  END IF
88  END IF
89 
90  IF (inputs%if_level_set) THEN
91  DO i = 1, m_max_c
92  DO k = 1, 2
93  DO int_nb = 1, inputs%nb_fluid - 1
94  level_setn_m1(int_nb,:,k,i) = level_setn(int_nb,:,k,i) &
95  - level_set_exact(int_nb,k,pp_mesh%rr,list_mode(i),time)
96  END DO
97  END DO
98  END DO
99  norm_err(3) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, level_setn_m1(1,:,:,:))
100  IF (rank==0) THEN
101  WRITE(10,*) 'Level set field#####################'
102  WRITE(10,*) 'L2 error on level set', norm_err(3)
103  WRITE(10,*) 'Level set field#####################'
104  END IF
105  END IF
106 
107  IF (inputs%numero_du_test_debug==24) THEN
108  DO i = 1, m_max_c
109  DO k= 1, 2
110  DO n = 1, SIZE(pp_mesh%rr,2)
111  IF (pp_mesh%rr(1,n).GT..5d0) THEN
112  pn_m1(n,k,i) = pn_m1(n,k,i)
113  ELSE
114  pn_m1(n,k,i) = 0.d0
115  END IF
116  END DO
117  END DO
118  END DO
119 
120  norm_err(4) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, pn_m1)
121  IF (rank==0) THEN
122  WRITE(10,*) 'L2 error on pressure outter obstacle = ', norm_err(4)
123  END IF
124  END IF
125 
126 
127  CASE(3:7,10,12)
128  DO i = 1, m_max_c
129  DO k =1, 6
130  hn1(:,k,i) = hn(:,k,i) - hexact(h_mesh, k, h_mesh%rr, list_mode(i), mu_h_field, time)
131  END DO
132  END DO
133 
134  IF (inputs%nb_dom_phi/=0) THEN
135  DO i = 1, m_max_c
136  DO k =1, 2
137  phin1(:,k,i) = phin(:,k,i) - phiexact(k, phi_mesh%rr, list_mode(i), inputs%mu_phi, time)
138  END DO
139  END DO
140  END IF
141 
142  norm = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn)
143  err = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn1)
144  norm_err(1) = err/norm
145 
146  norm = norm_sf(comm_one_d, 'sH1', h_mesh, list_mode, hn)
147  err = norm_sf(comm_one_d, 'curl', h_mesh, list_mode, hn1)
148  norm_err(2) = err/norm
149 
150  norm = norm_sf(comm_one_d, 'sH1', h_mesh, list_mode, bn)
151  err = norm_sf(comm_one_d, 'div', h_mesh, list_mode, bn)
152 
153  norm_err(3) = err/norm
154 
155  IF (inputs%nb_dom_phi/=0) THEN
156  norm = norm_sf(comm_one_d, 'sH1', phi_mesh, list_mode, phin)
157  err = norm_sf(comm_one_d, 'sH1', phi_mesh, list_mode, phin1)
158  norm_err(4) = err/norm
159  ELSE
160  norm_err(4) = inputs%norm_ref(4)
161  END IF
162 
163  IF (rank==0) THEN
164  WRITE(rank+10,*) 'Magnetic field #####################'
165  WRITE(rank+10,*) 'L2 norm of error on Hn = ', norm_err(1)
166  WRITE(rank+10,*) 'L2 norm of error on Curl(Hn) = ', norm_err(2)
167  WRITE(rank+10,*) 'L2 norm of Div(mu Hn) = ', norm_err(3)
168  IF (inputs%nb_dom_phi/=0) THEN
169  WRITE(10,*) 'Scal potential #####################'
170  WRITE(10,*) 'H1 norm of error on phin = ', norm_err(4)
171  END IF
172  END IF
173 
174  CASE(11)
175  norm_err(1) = norm_sf(comm_one_d_ns, 'sH1', vv_mesh, list_mode, un)
176  norm_err(2) = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn)
177  norm_err(3) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, pn)
178  norm_err(4) = norm_sf(comm_one_d_ns, 'L2', temp_mesh, list_mode, tempn)
179  IF (rank==0) THEN
180  WRITE(10,*) '######################################'
181  WRITE(10,*) 'H1 norm on velocity = ', norm_err(1)
182  WRITE(10,*) 'L2 norm of Hn = ', norm_err(2)
183  WRITE(10,*) 'L2 norm of pressure = ', norm_err(3)
184  WRITE(10,*) 'L2 norm of temperature = ', norm_err(4)
185  WRITE(10,*) '######################################'
186  END IF
187 
188  CASE(13)
189  norm_err(1) = norm_sf(comm_one_d_ns, 'sH1', vv_mesh, list_mode, un)
190  norm_err(2) = norm_sf(comm_one_d, 'div', h_mesh, list_mode, hn)
191  norm_err(3) = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn)
192  norm_err(4) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, pn)
193  IF (rank==0) THEN
194  WRITE(10,*) '######################################'
195  WRITE(10,*) 'H1 norm on velocity = ', norm_err(1)
196  WRITE(10,*) 'L2 norm of div(Hn) = ', norm_err(2)
197  WRITE(10,*) 'L2 norm of Hn = ', norm_err(3)
198  WRITE(10,*) 'L2 norm of pressure = ', norm_err(4)
199  WRITE(10,*) '######################################'
200  END IF
201 
202  CASE(14)
203  IF (rank==0) THEN
204  READ(10,*) norm_err(1)
205  READ(10,*) norm_err(2)
206  READ(11,*) norm_err(3)
207  READ(11,*) norm_err(4)
208  END IF
209 
210  CASE(15)
211 
212  norm_err_loc = 0.d0
213  DO i = 1, m_max_c
214  norm_err_loc(rank+1)= 0.5*(norme_l2_champ_par(comm_one_d_ns(1), &
215  vv_mesh, list_mode(i:i), un(:,:,i:i)))**2
216  END DO
217  norm_err_loc(4)= norm_sf(comm_one_d_ns, 'div', vv_mesh, list_mode, un)/&
218  norm_sf(comm_one_d_ns, 'sH1', vv_mesh, list_mode, un)
219  CALL mpi_allreduce(norm_err_loc,norm_err,4,mpi_double_precision, &
220  mpi_max, comm_one_d(2), code)
221 
222  IF (rank==0) THEN
223  WRITE(10,*) '######################################'
224  WRITE(10,*) 'L2 norm of mode 0 = ', norm_err(1)
225  WRITE(10,*) 'L2 norm of mode 1 = ', norm_err(2)
226  WRITE(10,*) 'L2 norm of mode 2 = ', norm_err(3)
227  WRITE(10,*) 'Relative L2 norm of div.= ', norm_err(4)
228  WRITE(10,*) '######################################'
229  END IF
230 
231  CASE(16)
232  norm_err(1) = 0.5*norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un)**2
233  CALL angular_momentum(vv_mesh, list_mode, un, norm_err(2:4))
234  IF (rank==0) THEN
235  WRITE(10,*) '######################################'
236  WRITE(10,*) 'Total kinetic energy = ', norm_err(1)
237  WRITE(10,*) 'Mx = ', norm_err(2)
238  WRITE(10,*) 'My = ', norm_err(3)
239  WRITE(10,*) 'Mz = ', norm_err(4)
240  WRITE(10,*) '######################################'
241  END IF
242 
243  CASE(17,22,27)
244  norm = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, bn)
245  err = norm_sf(comm_one_d, 'div', h_mesh, list_mode, bn)
246  norm_err(1) = err/norm
247  DO i = 1, m_max_c
248  DO k = 1, 6
249  hn1(:,k,i) = hexact(h_mesh, k, h_mesh%rr, list_mode(i), mu_h_field, time)
250  hn1(:,k,i) = hn1(:,k,i) - hn(:,k,i)
251  ENDDO
252  ENDDO
253  err = norm_sf(comm_one_d, 'curl', h_mesh, list_mode, hn1)
254  norm = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn)
255  norm_err(2) = err/norm
256  err = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn1)
257  norm_err(3) = err/norm
258  DO i = 1, m_max_c
259  DO k = 1, 2
260  phin1(:,k,i) = phin(:,k,i) - phiexact(k, phi_mesh%rr, list_mode(i), inputs%mu_phi, time)
261  END DO
262  END DO
263  err = norm_sf(comm_one_d, 'L2', phi_mesh, list_mode, phin1)
264  norm = norm_sf(comm_one_d, 'L2', phi_mesh, list_mode, phin)
265  norm_err(4) = err/norm
266  IF (rank==0) THEN
267  WRITE(10,*) '########################################################'
268  WRITE(10,*) ' L2-norm on div of B err/norm =', norm_err(1)
269  WRITE(10,*) ' L2-norm on curl of (H-Hexact) err/norm =', norm_err(2)
270  WRITE(10,*) ' L2-norm of (H-Hexact) err/norm =', norm_err(3)
271  WRITE(10,*) ' L2-norm on phi err/norm =', norm_err(4)
272  WRITE(10,*) '########################################################'
273  END IF
274 
275  CASE(18,23,26)
276  norm = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, bn)
277  err = norm_sf(comm_one_d, 'div', h_mesh, list_mode, bn)
278  norm_err(1) = err/norm
279  DO i = 1, m_max_c
280  DO k = 1, 6
281  hn1(:,k,i) = hexact(h_mesh, k, h_mesh%rr, list_mode(i), mu_h_field, time)
282  hn1(:,k,i) = hn1(:,k,i) - hn(:,k,i)
283  ENDDO
284  ENDDO
285  err = norm_sf(comm_one_d, 'curl', h_mesh, list_mode, hn1)
286  norm = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn)
287  norm_err(2) = err/norm
288  err = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn1)
289  norm_err(3) = err/norm
290  norm_err(4) = norm
291  IF (rank==0) THEN
292  WRITE(10,*) '########################################################'
293  WRITE(10,*) ' L2-norm on div of B err/norm =', norm_err(1)
294  WRITE(10,*) ' L2-norm on curl of (H-Hexact) err/norm =', norm_err(2)
295  WRITE(10,*) ' L2-norm of (H-Hexact) err/norm =', norm_err(3)
296  WRITE(10,*) ' L2-norm of H =', norm_err(4)
297  WRITE(10,*) '########################################################'
298  END IF
299 
300  CASE(21)
301  DO i = 1, m_max_c
302  DO k= 1, 6
303  un_m1(:,k,i) = un(:,k,i) - vv_exact(k,vv_mesh%rr,list_mode(i),time)
304  hn1(:,k,i) = hn(:,k,i) - hexact(h_mesh, k, h_mesh%rr, list_mode(i), mu_h_field, time)
305  END DO
306  DO k= 1, 2
307  pn_m1(:,k,i) = pn(:,k,i) - pp_exact(k,pp_mesh%rr,list_mode(i),time)
308  DO int_nb = 1, inputs%nb_fluid - 1
309  level_setn_m1(int_nb,:,k,i) = level_setn(int_nb,:,k,i) &
310  - level_set_exact(int_nb,k,pp_mesh%rr,list_mode(i),time)
311  END DO
312  END DO
313  IF (list_mode(i) == 0) THEN
314  CALL moy(comm_one_d(1),pp_mesh, pn_m1(:,1,i),moyenne)
315  pn_m1(:,1,i) = pn_m1(:,1,i) - moyenne
316  ENDIF
317  END DO
318  norm_err(1) = norm_sf(comm_one_d_ns, 'L2', vv_mesh, list_mode, un_m1)
319  norm_err(2) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, pn_m1)
320  norm_err(3) = norm_sf(comm_one_d_ns, 'L2', pp_mesh, list_mode, level_setn_m1(1,:,:,:))
321  err = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn1)
322  norm_err(4) = err
323  IF (rank==0) THEN
324  WRITE(10,*) 'Velocity field #####################'
325  WRITE(10,*) 'L2 error on velocity = ', norm_err(1)
326  WRITE(10,*) 'Pressure field #####################'
327  WRITE(10,*) 'L2 error on pressure = ', norm_err(2)
328  WRITE(10,*) 'Level Set #####################'
329  WRITE(10,*) 'L2 error on level set = ', norm_err(3)
330  WRITE(10,*) 'Magnetic field #####################'
331  WRITE(10,*) 'L2 error on Hn = ', norm_err(4)
332  END IF
333 
334 
335  CASE(28)
336  err = norm_sf(comm_one_d, 'div', vv_mesh, list_mode, un)
337  norm = norm_sf(comm_one_d, 'H1', vv_mesh, list_mode, un)
338  norm_err(1) = err
339  norm_err(2) = err/norm
340  norm = norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un)
341  norm_err(3) = norm
342  norm = norm_sf(comm_one_d, 'H1', pp_mesh, list_mode, pn)
343  norm_err(4) = norm
344 
345  IF (rank==0) THEN
346  WRITE(10,*) '########################################################'
347  WRITE(10,*) ' L2-norm on div of u =', norm_err(1)
348  WRITE(10,*) ' L2-norm of div of u err/norm =', norm_err(2)
349  WRITE(10,*) ' L2-norm of u =', norm_err(3)
350  WRITE(10,*) ' H1-norm of p =', norm_err(4)
351  WRITE(10,*) '########################################################'
352  END IF
353 
354  CASE(29)
355 
356  err = norm_sf(comm_one_d, 'div', vv_mesh, list_mode, un)
357  norm = norm_sf(comm_one_d, 'H1', vv_mesh, list_mode, un)
358  norm_err(1) = err/norm
359  norm = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, bn)
360  err = norm_sf(comm_one_d, 'div', h_mesh, list_mode, bn)
361  norm_err(2) = err/norm
362  norm_err(3) = norm
363  err= dot_product_sf(comm_one_d, h_mesh, list_mode, hn, bn)
364  norm_err(4) = 0.5*err
365  IF (rank==0) THEN
366  WRITE(10,*) '########################################################'
367  WRITE(10,*) ' L2-norm of div of u err/norm =', norm_err(1)
368  WRITE(10,*) ' L2-norm on div of B err/norm =', norm_err(2)
369  WRITE(10,*) ' L2-norm of of B =', norm_err(3)
370  WRITE(10,*) ' integral of 0.5*B.H =', norm_err(4)
371  WRITE(10,*) '########################################################'
372  END IF
373 
374  CASE(30)
375 
376  DO type=1,6
377  DO i=1,size(list_mode)
378  un_ex(:,type,i) = vv_exact(type,vv_mesh%rr,list_mode(i),time)
379  END DO
380  END DO
381  un_error = un - un_ex
382  norm_err(1) = norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_error) / norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_ex)
383  DO type=1,2
384  DO i=1,size(list_mode)
385  pn_ex(:,type,i) = pp_exact(type,pp_mesh%rr,list_mode(i),time)
386  END DO
387  END DO
388  pn_error = pn - pn_ex
389  norm_err(2) = norm_sf(comm_one_d, 'L2', pp_mesh, list_mode, pn_error)
390  DO type=1,2
391  DO i=1,size(list_mode)
392  tempn_ex(:,type,i) = temperature_exact(type,temp_mesh%rr,list_mode(i),time)
393  END DO
394  END DO
395  tempn_error = tempn - tempn_ex
396  norm_err(3) = norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_error) / &
397  norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_ex)
398  norm_err(4) = norm_sf(comm_one_d, 'H1', temp_mesh, list_mode, tempn_error) / &
399  norm_sf(comm_one_d, 'H1', temp_mesh, list_mode, tempn_ex)
400  IF (rank==0) THEN
401  WRITE(10,*) '########################################################'
402  WRITE(10,*) ' L2-norm of error on u / L2-norm of u exact =', norm_err(1)
403  WRITE(10,*) ' L2-norm of error on p =', norm_err(2)
404  WRITE(10,*) ' L2-norm of error on T / L2-norm of T exact =', norm_err(3)
405  WRITE(10,*) ' H1-norm of error on T / H1-norm of T exact =', norm_err(4)
406  WRITE(10,*) '########################################################'
407  END IF
408 
409  CASE(31,32)
410 
411  DO type=1,6
412  DO i=1,size(list_mode)
413  un_ex(:,type,i) = vv_exact(type,vv_mesh%rr,list_mode(i),time)
414  END DO
415  END DO
416  un_error = un - un_ex
417  norm_err(1) = norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_error) / norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_ex)
418  DO type=1,2
419  DO i=1,size(list_mode)
420  pn_ex(:,type,i) = pp_exact(type,pp_mesh%rr,list_mode(i),time)
421  END DO
422  END DO
423  pn_error = pn - pn_ex
424  norm_err(2) = norm_sf(comm_one_d, 'L2', pp_mesh, list_mode, pn_error) / norm_sf(comm_one_d, 'L2', pp_mesh, list_mode, pn_ex)
425  DO type=1,2
426  DO i=1,size(list_mode)
427  tempn_ex(:,type,i) = temperature_exact(type,temp_mesh%rr,list_mode(i),time)
428  END DO
429  END DO
430  tempn_error = tempn - tempn_ex
431  norm_err(3) = norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_error) / &
432  norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_ex)
433  norm_err(4) = norm_sf(comm_one_d, 'H1', temp_mesh, list_mode, tempn_error) / &
434  norm_sf(comm_one_d, 'H1', temp_mesh, list_mode, tempn_ex)
435  IF (rank==0) THEN
436  WRITE(10,*) '########################################################'
437  WRITE(10,*) ' L2-norm of error on u / L2-norm of u exact =', norm_err(1)
438  WRITE(10,*) ' L2-norm of error on p / L2-norm of p exact =', norm_err(2)
439  WRITE(10,*) ' L2-norm of error on T / L2-norm of T exact =', norm_err(3)
440  WRITE(10,*) ' H1-norm of error on T / H1-norm of T exact =', norm_err(4)
441  WRITE(10,*) '########################################################'
442  END IF
443 
444  CASE(33)
445 
446  DO type=1,6
447  DO i=1,size(list_mode)
448  un_ex(:,type,i) = vv_exact(type,vv_mesh%rr,list_mode(i),time)
449  END DO
450  END DO
451  un_error = un - un_ex
452  norm_err(1) = norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_error) / norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_ex)
453  DO type=1,2
454  DO i=1,size(list_mode)
455  pn_ex(:,type,i) = pp_exact(type,pp_mesh%rr,list_mode(i),time)
456  END DO
457  END DO
458  pn_error = pn - pn_ex
459  norm_err(2) = norm_sf(comm_one_d, 'L2', pp_mesh, list_mode, pn_error)
460  DO type=1,2
461  DO i=1,size(list_mode)
462  tempn_ex(:,type,i) = temperature_exact(type,temp_mesh%rr,list_mode(i),time)
463  END DO
464  END DO
465  tempn_error = tempn - tempn_ex
466  norm_err(3) = norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_error) / &
467  norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_ex)
468  DO type=1,6
469  DO i=1,size(list_mode)
470  hn_ex(:,type,i) = hexact(h_mesh,type,h_mesh%rr,list_mode(i),mu_h_field,time)
471  END DO
472  END DO
473  hn_error = hn - hn_ex
474  norm_err(4) = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn_error) / norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn_ex)
475  IF (rank==0) THEN
476  WRITE(10,*) '########################################################'
477  WRITE(10,*) ' L2-norm of error on u / L2-norm of u exact =', norm_err(1)
478  WRITE(10,*) ' L2-norm of error on p =', norm_err(2)
479  WRITE(10,*) ' L2-norm of error on T / L2-norm of T exact =', norm_err(3)
480  WRITE(10,*) ' L2-norm of error on H / L2-norm of H exact =', norm_err(4)
481  WRITE(10,*) '########################################################'
482  END IF
483 
484  CASE(34)
485 
486  DO type=1,6
487  DO i=1,size(list_mode)
488  un_ex(:,type,i) = vv_exact(type,vv_mesh%rr,list_mode(i),time)
489  END DO
490  END DO
491  un_error = un - un_ex
492  norm_err(1) = norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_error) / norm_sf(comm_one_d, 'L2', vv_mesh, list_mode, un_ex)
493  DO type=1,2
494  DO i=1,size(list_mode)
495  pn_ex(:,type,i) = pp_exact(type,pp_mesh%rr,list_mode(i),time)
496  END DO
497  END DO
498  pn_error = pn - pn_ex
499  norm_err(2) = norm_sf(comm_one_d, 'L2', pp_mesh, list_mode, pn_error) / norm_sf(comm_one_d, 'L2', pp_mesh, list_mode, pn_ex)
500  DO type=1,2
501  DO i=1,size(list_mode)
502  tempn_ex(:,type,i) = temperature_exact(type,temp_mesh%rr,list_mode(i),time)
503  END DO
504  END DO
505  tempn_error = tempn - tempn_ex
506  norm_err(3) = norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_error) / &
507  norm_sf(comm_one_d, 'L2', temp_mesh, list_mode, tempn_ex)
508  DO type=1,6
509  DO i=1,size(list_mode)
510  hn_ex(:,type,i) = hexact(h_mesh,type,h_mesh%rr,list_mode(i),mu_h_field,time)
511  END DO
512  END DO
513  hn_error = hn - hn_ex
514  norm_err(4) = norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn_error) / norm_sf(comm_one_d, 'L2', h_mesh, list_mode, hn_ex)
515  IF (rank==0) THEN
516  WRITE(10,*) '########################################################'
517  WRITE(10,*) ' L2-norm of error on u / L2-norm of u exact =', norm_err(1)
518  WRITE(10,*) ' L2-norm of error on p / L2-norm of p exact =', norm_err(2)
519  WRITE(10,*) ' L2-norm of error on T / L2-norm of T exact =', norm_err(3)
520  WRITE(10,*) ' L2-norm of error on H / L2-norm of H exact =', norm_err(4)
521  WRITE(10,*) '########################################################'
522  END IF
523 
524  CASE default
525  CALL error_petsc(' BUG in post_proc_test: We should not be here')
526 
527  END SELECT
528 
529  IF (rank==0) THEN
530  ! Compare with references
531 !!$ IF (MAXVAL(ABS(inputs%norm_ref-norm_err)/inputs%norm_ref)<error_on_tests) THEN
532  IF (maxval(abs(inputs%norm_ref-norm_err))<1.d-8) THEN
533  WRITE(57,'(A,I2,A)') 'test #',inputs%numero_du_test_debug,' OK'
534  ELSE
535 !!$ WRITE(57,'(A,I2,2x,e15.7)') 'Problem with test #', inputs%numero_du_test_debug, &
536 !!$ MAXVAL(ABS(inputs%norm_ref-norm_err)/inputs%norm_ref)
537  WRITE(57,'(A,I2,2x,e15.7)') 'Problem with test #', inputs%numero_du_test_debug, &
538  maxval(abs(inputs%norm_ref-norm_err))
539  END IF
540  ! End compare
541  END IF
542 
543  END SUBROUTINE post_proc_test
544  !---------------------------------------------------------------------------
545 
546 END MODULE post_processing_debug
subroutine angular_momentum(mesh, list_mode, field, moments_out)
Definition: tn_axi.f90:247
Definition: tn_axi.f90:5
real(kind=8) function, dimension(size(rr, 2)), public pp_exact(TYPE, rr, m, t)
subroutine, public post_proc_test(vv_mesh, pp_mesh, temp_mesh, H_mesh, phi_mesh, list_mode, un, pn, Hn, Bn, phin, tempn, level_setn, mu_H_field, time, m_max_c, comm_one_d, comm_one_d_ns)
subroutine moy(mesh, p, RESULT)
real(kind=8) function, dimension(size(rr, 2)), public temperature_exact(TYPE, rr, m, t)
real(kind=8) function, dimension(size(rr, 2)), public vv_exact(TYPE, rr, m, t)
real(kind=8) function, dimension(size(rr, 2)), public phiexact(TYPE, rr, m, mu_phi, t)
real(kind=8) function norm_sf(communicator, norm_type, mesh, list_mode, v)
Definition: tn_axi.f90:38
real(kind=8) function norme_l2_champ_par(communicator, mesh, list_mode, v)
Definition: tn_axi.f90:179
subroutine error_petsc(string)
Definition: my_util.f90:15
real(kind=8) function dot_product_sf(communicator, mesh, list_mode, v, w)
Definition: tn_axi.f90:9
real(kind=8) function, dimension(size(rr, 2)), public level_set_exact(interface_nb, TYPE, rr, m, t)
real(kind=8) function, dimension(size(rr, 2)), public hexact(H_mesh, TYPE, rr, m, mu_H_field, t)