Choral
quad_mod.f90
Go to the documentation of this file.
1 
30 
31 module quad_mod
32 
34  use choral_env
35  use real_type
36  use io, only: read
37  use cell_mod
38 
39  implicit none
40  private
41 
42 
43  ! %----------------------------------------%
44  ! | |
45  ! | PUBLIC DATA |
46  ! | |
47  ! %----------------------------------------%
48 
49  public :: quad_init ! TESTED
50 
51  public :: quad_name
52  public :: quad_nbnodes
53  public :: quad_geo
54  public :: quad_order
55  public :: quad_dim
56  public :: quad_coord
57  public :: quad_wgt
58 
59  ! %----------------------------------------%
60  ! | |
61  ! | QUAD ARRAYS |
62  ! | |
63  ! %----------------------------------------%
64 
68 
70  character(LEN=13), dimension(0:QUAD_TOT_NB) :: quad_name
71 
73  integer , dimension(QUAD_TOT_NB):: quad_nbnodes
74 
76  integer , dimension(QUAD_TOT_NB):: quad_geo
77 
79  integer , dimension(QUAD_TOT_NB):: quad_dim
80 
82  integer , dimension(QUAD_TOT_NB):: quad_order
83 
85  type(r_2d), dimension(QUAD_TOT_NB):: quad_coord
86 
88  type(r_1d), dimension(QUAD_TOT_NB):: quad_wgt
89 
90 
91 contains
92 
93 
96  subroutine quad_init(b)
97  logical, intent(in) :: b
98 
99  integer :: qt, geo, dim, nn
100 
101  if (b) write(*,*) "quad_mod : quad_init"
102 
103  !! set quad name, geo, order, nbNodes
104  !!
105  call def_quad_xxx()
106 
107  !! set quad dimension
108  !!
109  do qt=1, quad_tot_nb
110  geo = quad_geo(qt)
111  dim = cell_dim(geo)
112  quad_dim(qt) = dim
113  end do
114 
115  !! allocate memory to store node coordinates
116  !! and weights
117  !!
118  do qt=1, quad_tot_nb
119 
120  nn = quad_nbnodes(qt)
121  dim = quad_dim(qt)
122 
123  quad_coord(qt) = r_2d( dim, nn)
124  quad_wgt(qt) = r_1d( nn)
125 
126  end do
127 
128  !! set GUAD_COORD AND QUAD_WGT
129  !!
130  call def_quad_data()
131 
132  end subroutine quad_init
133 
134  subroutine def_quad_xxx()
136  quad_name(quad_none ) = "VOID"
137  quad_name(quad_gauss_edg_1 ) = "Gauss EDG 1 "
138  quad_name(quad_gauss_edg_2 ) = "Gauss EDG 2 "
139  quad_name(quad_gauss_edg_3 ) = "Gauss EDG 3 "
140  quad_name(quad_gauss_edg_4 ) = "Gauss EDG 4 "
141  quad_name(quad_gauss_trg_1 ) = "Gauss TRG 1 "
142  quad_name(quad_gauss_trg_3 ) = "Gauss TRG 3 "
143  quad_name(quad_gauss_trg_6 ) = "Gauss TRG 6 " ! STRANG_5
144  quad_name(quad_gauss_trg_12) = "Gauss TRG 12 " ! STRANG_9
145  quad_name(quad_gauss_trg_13) = "Gauss TRG 13 " ! STRANG_10
146  quad_name(quad_gauss_trg_19) = "Gauss TRG 19 " ! TOMS584_19
147  quad_name(quad_gauss_trg_28) = "Gauss TRG 28 " ! TOMS612_28
148  quad_name(quad_gauss_trg_37) = "Gauss TRG 37 " ! TOMS706_37
149  quad_name(quad_gauss_tet_1 ) = "Gauss TET 1 "
150  quad_name(quad_gauss_tet_4 ) = "Gauss TET 4 " ! keast1 1
151  quad_name(quad_gauss_tet_15) = "Gauss TET 15 " ! keast1 6
152  quad_name(quad_gauss_tet_31) = "Gauss TET 31 " ! keast1 8
153  quad_name(quad_gauss_tet_45) = "Gauss TET 45 " ! keast1 9
154 
172 
190 
201  quad_order(quad_gauss_trg_28) = 10 ! under converging /= 11
202  quad_order(quad_gauss_trg_37) = 10 ! under converging /= 13
208 
209  end subroutine def_quad_xxx
210 
211 
212  subroutine def_quad_data()
214  real(RP) :: a
215  integer :: qt
216 
217  qt = quad_gauss_edg_1
218  quad_wgt(qt)%y = 1.0_rp
219  quad_coord(qt)%y = 0.5_rp
220 
221  qt = quad_gauss_edg_2
222  quad_wgt(qt)%y = 0.5_rp
223  quad_coord(qt)%y(1,1) = 0.21132486540518702_rp
224  quad_coord(qt)%y(1,2) = 0.788675134594813_rp
225 
226  qt = quad_gauss_edg_3
227  quad_wgt(qt)%y(1) = 0.44444444444444444444444444444444444_rp
228  quad_wgt(qt)%y(2) = 0.27777777777777777777777777777777778_rp
229  quad_wgt(qt)%y(3) = 0.27777777777777777777777777777777778_rp
230  quad_coord(qt)%y(1,1) = 0.5_rp
231  quad_coord(qt)%y(1,2) = 0.11270166537925852_rp
232  quad_coord(qt)%y(1,3) = 0.8872983346207415_rp
233 
234  qt = quad_gauss_edg_4
235  quad_wgt(qt)%y(1) = 0.17392742256872692868653197461099970_rp
236  quad_wgt(qt)%y(2) = 0.17392742256872692868653197461099970_rp
237  quad_wgt(qt)%y(3) = 0.32607257743127307131346802538900030_rp
238  quad_wgt(qt)%y(4) = 0.32607257743127307131346802538900030_rp
239  a = 0.06943184420297371238802675555359525_rp
240  quad_coord(qt)%y(1,1) = a
241  quad_coord(qt)%y(1,2) = 1._rp - a
242  a = 0.33000947820757186759866712044837766_rp
243  quad_coord(qt)%y(1,3) = a
244  quad_coord(qt)%y(1,4) = 1._rp-a
245 
246  qt = quad_gauss_trg_1
247  quad_wgt(qt)%y = 0.5_rp
248  quad_coord(qt)%y(:,1) = re(1, 3)
249 
250  qt = quad_gauss_trg_3
251  quad_wgt(qt)%y = re(1, 6)
252  quad_coord(qt)%y(:,1) = re(1, 6)
253  quad_coord(qt)%y(:,2) = (/ re(2, 3), re(1, 6)/)
254  quad_coord(qt)%y(:,3) = (/ re(1, 6), re(2, 3)/)
255 
256  qt = quad_gauss_trg_6
257  call read(quad_wgt(qt)%y, &
258  & trim(choral_dir)//'SRC/quad_data/strang5_w.txt')
259  quad_wgt(qt)%y = quad_wgt(qt)%y * 0.5_rp
260  call read(quad_coord(qt)%y, &
261  & trim(choral_dir)//'SRC/quad_data/strang5_x.txt')
262 
263 
264  qt = quad_gauss_trg_12
265  call read(quad_wgt(qt)%y, &
266  & trim(choral_dir)//'SRC/quad_data/strang9_w.txt')
267  quad_wgt(qt)%y = quad_wgt(qt)%y*0.5_rp
268  call read(quad_coord(qt)%y, &
269  & trim(choral_dir)//'SRC/quad_data/strang9_x.txt')
270 
271  qt = quad_gauss_trg_13
272  call read(quad_wgt(qt)%y, &
273  & trim(choral_dir)//'SRC/quad_data/strang10_w.txt')
274  quad_wgt(qt)%y = quad_wgt(qt)%y*0.5_rp
275  call read(quad_coord(qt)%y, &
276  & trim(choral_dir)//'SRC/quad_data/strang10_x.txt')
277 
278  qt = quad_gauss_trg_19
279  call read(quad_wgt(qt)%y, &
280  & trim(choral_dir)//'SRC/quad_data/TOMS584_19_w.txt')
281  quad_wgt(qt)%y = quad_wgt(qt)%y*0.5_rp
282  call read(quad_coord(qt)%y, &
283  & trim(choral_dir)//'SRC/quad_data/TOMS584_19_x.txt')
284 
285 
286  qt = quad_gauss_trg_28
287  call read(quad_wgt(qt)%y, &
288  & trim(choral_dir)//'SRC/quad_data/toms612_28_w.txt')
289  quad_wgt(qt)%y = quad_wgt(qt)%y*0.5_rp
290  call read(quad_coord(qt)%y, &
291  & trim(choral_dir)//'SRC/quad_data/toms612_28_x.txt')
292 
293  qt = quad_gauss_trg_37
294  call read(quad_wgt(qt)%y, &
295  & trim(choral_dir)//'SRC/quad_data/toms706_37_w.txt')
296  quad_wgt(qt)%y = quad_wgt(qt)%y*0.5_rp
297  call read(quad_coord(qt)%y, &
298  & trim(choral_dir)//'SRC/quad_data/toms706_37_x.txt')
299 
300 
301  qt = quad_gauss_tet_1
302  quad_wgt(qt)%y = re(1, 6)
303  quad_coord(qt)%y = re(1, 4)
304 
305  qt = quad_gauss_tet_4
306  call read(quad_wgt(qt)%y, &
307  & trim(choral_dir)//'SRC/quad_data/keast1_w.txt')
308  quad_wgt(qt)%y = quad_wgt(qt)%y / 6._rp
309  call read(quad_coord(qt)%y, &
310  & trim(choral_dir)//'SRC/quad_data/keast1_x.txt')
311 
312  qt = quad_gauss_tet_15
313  call read(quad_wgt(qt)%y, &
314  & trim(choral_dir)//'SRC/quad_data/keast6_w.txt')
315  quad_wgt(qt)%y = quad_wgt(qt)%y / 6._rp
316  call read(quad_coord(qt)%y, &
317  & trim(choral_dir)//'SRC/quad_data/keast6_x.txt')
318 
319  qt = quad_gauss_tet_31
320  call read(quad_wgt(qt)%y, &
321  & trim(choral_dir)//'SRC/quad_data/keast8_w.txt')
322  quad_wgt(qt)%y = quad_wgt(qt)%y / 6._rp
323  call read(quad_coord(qt)%y, &
324  & trim(choral_dir)//'SRC/quad_data/keast8_x.txt')
325 
326  qt = quad_gauss_tet_45
327  call read(quad_wgt(qt)%y, &
328  & trim(choral_dir)//'SRC/quad_data/keast9_w.txt')
329  quad_wgt(qt)%y = quad_wgt(qt)%y / 6._rp
330  call read(quad_coord(qt)%y, &
331  & trim(choral_dir)//'SRC/quad_data/keast9_x.txt')
332 
333  end subroutine def_quad_data
334 
335 
336 end module quad_mod
integer, parameter quad_gauss_trg_28
integer, parameter cell_tet
Tetrahedron.
integer, parameter quad_gauss_tet_15
type(r_1d), dimension(quad_tot_nb), public quad_wgt
quad weights
Definition: quad_mod.f90:88
integer, parameter quad_gauss_edg_4
QUADRATURE RULES ON REFERENCE CELLS
Definition: quad_mod.f90:31
subroutine def_quad_data()
Definition: quad_mod.f90:213
integer, parameter quad_gauss_tet_4
integer, parameter quad_gauss_trg_3
integer, parameter quad_none
subroutine, public quad_init(b)
initialise QUAD_XXX arrays:
Definition: quad_mod.f90:97
integer, parameter cell_edg
Edge (line segment)
integer, parameter quad_gauss_tet_1
integer, dimension(cell_tot_nb), public cell_dim
Dimension for each cell type.
Definition: cell_mod.f90:130
integer, parameter quad_gauss_tet_45
integer, dimension(quad_tot_nb), public quad_dim
Reference cell dimension for each quad method.
Definition: quad_mod.f90:79
conversion integers or rational to real
Definition: real_type.F90:153
integer, dimension(quad_tot_nb), public quad_order
Order for each quad method.
Definition: quad_mod.f90:82
character(len=100), parameter, public choral_dir
CHORAL_DIR = path to choral.
Definition: choral_env.f90:27
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
Definition: real_type.F90:33
IO: module for input/output
Definition: io.f90:8
type(r_2d), dimension(quad_tot_nb), public quad_coord
quad node coordinates
Definition: quad_mod.f90:85
integer, parameter quad_gauss_edg_3
integer, parameter quad_tot_nb
Total number of quadrature rules.
integer, parameter cell_trg
Triangle.
integer, dimension(quad_tot_nb), public quad_nbnodes
Number of nodes for each quad method.
Definition: quad_mod.f90:73
CHORAL CONSTANTS
integer, parameter quad_gauss_trg_1
subroutine def_quad_xxx()
Definition: quad_mod.f90:135
R_1D: this type will allow to define arrays of 1D real arrays.
Definition: real_type.F90:106
integer, parameter quad_gauss_trg_19
integer, parameter quad_gauss_trg_37
DEFINITION OF DIRECTORY PATHS
Definition: choral_env.f90:16
integer, parameter quad_gauss_edg_2
integer, parameter quad_gauss_trg_13
read from file
Definition: io.f90:39
integer, parameter quad_gauss_edg_1
integer, parameter quad_gauss_tet_31
character(len=13), dimension(0:quad_tot_nb), public quad_name
QUAD_ARRAY.
Definition: quad_mod.f90:70
R_2D: this type will allow to define arrays of 2D real arrays.
Definition: real_type.F90:114
integer, parameter quad_gauss_trg_6
integer, parameter quad_gauss_trg_12
DEFINITION OF GEOMETRICAL CELLS (for meshes)
Definition: cell_mod.f90:82
integer, dimension(quad_tot_nb), public quad_geo
Reference cell geometry for each quad method.
Definition: quad_mod.f90:76