geometry_base.f90
Go to the documentation of this file.
1!#############################################################################
2!# #
3!# fosite - 3D hydrodynamical simulation program #
4!# module: geometry_generic.f90 #
5!# #
6!# Copyright (C) 2007-2018 #
7!# Tobias Illenseer <tillense@astrophysik.uni-kiel.de> #
8!# Jannes Klee <jklee@astrophysik.uni-kiel.de> #
9!# #
10!# This program is free software; you can redistribute it and/or modify #
11!# it under the terms of the GNU General Public License as published by #
12!# the Free Software Foundation; either version 2 of the License, or (at #
13!# your option) any later version. #
14!# #
15!# This program is distributed in the hope that it will be useful, but #
16!# WITHOUT ANY WARRANTY; without even the implied warranty of #
17!# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, GOOD TITLE or #
18!# NON INFRINGEMENT. See the GNU General Public License for more #
19!# details. #
20!# #
21!# You should have received a copy of the GNU General Public License #
22!# along with this program; if not, write to the Free Software #
23!# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. #
24!# #
25!#############################################################################
26
27!----------------------------------------------------------------------------!
35!----------------------------------------------------------------------------!
40 USE common_dict
41 !--------------------------------------------------------------------------!
42 PRIVATE
43 REAL, PARAMETER :: pi = 3.1415926535897932384626433832795028842
44! INTERFACE GetScale
45! MODULE PROCEDURE GetScale1, GetScale2
46! END INTERFACE
47! INTERFACE SetScale
48! MODULE PROCEDURE SetScale1, SetScale2, SetScale3
49! END INTERFACE
50
51 TYPE, ABSTRACT, EXTENDS (logging_base) :: geometry_base
52 REAL,DIMENSION(3) :: geoparam
53 INTEGER, PRIVATE :: azimuthindex = 0
54 CONTAINS
55
56! PRIVATE
57 PROCEDURE :: scalefactors_0
58 procedure(scalefactors_1), DEFERRED :: scalefactors_1
59 procedure(scalefactors_2), DEFERRED :: scalefactors_2
60 procedure(scalefactors_3), DEFERRED :: scalefactors_3
61 procedure(scalefactors_4), DEFERRED :: scalefactors_4
62 PROCEDURE :: radius_0
63 procedure(radius_1), DEFERRED :: radius_1
64 procedure(radius_2), DEFERRED :: radius_2
65 procedure(radius_3), DEFERRED :: radius_3
66 procedure(radius_4), DEFERRED :: radius_4
67 PROCEDURE :: positionvector_0
68 procedure(positionvector_1), DEFERRED :: positionvector_1
69 procedure(positionvector_2), DEFERRED :: positionvector_2
70 procedure(positionvector_3), DEFERRED :: positionvector_3
71 procedure(positionvector_4), DEFERRED :: positionvector_4
73 procedure(convert2cartesian_coords_1), DEFERRED :: convert2cartesian_coords_1
74 procedure(convert2cartesian_coords_2), DEFERRED :: convert2cartesian_coords_2
75 procedure(convert2cartesian_coords_3), DEFERRED :: convert2cartesian_coords_3
76 procedure(convert2cartesian_coords_4), DEFERRED :: convert2cartesian_coords_4
78 procedure(convert2curvilinear_coords_1), DEFERRED :: convert2curvilinear_coords_1
79 procedure(convert2curvilinear_coords_2), DEFERRED :: convert2curvilinear_coords_2
80 procedure(convert2curvilinear_coords_3), DEFERRED :: convert2curvilinear_coords_3
81 procedure(convert2curvilinear_coords_4), DEFERRED :: convert2curvilinear_coords_4
83 procedure(convert2cartesian_vectors_1), DEFERRED :: convert2cartesian_vectors_1
84 procedure(convert2cartesian_vectors_2), DEFERRED :: convert2cartesian_vectors_2
85 procedure(convert2cartesian_vectors_3), DEFERRED :: convert2cartesian_vectors_3
86 procedure(convert2cartesian_vectors_4), DEFERRED :: convert2cartesian_vectors_4
88 procedure(convert2curvilinear_vectors_1), DEFERRED :: convert2curvilinear_vectors_1
89 procedure(convert2curvilinear_vectors_2), DEFERRED :: convert2curvilinear_vectors_2
90 procedure(convert2curvilinear_vectors_3), DEFERRED :: convert2curvilinear_vectors_3
91 procedure(convert2curvilinear_vectors_4), DEFERRED :: convert2curvilinear_vectors_4
92 PROCEDURE :: setscale1
93 PROCEDURE :: setscale2
94 PROCEDURE :: setscale3
95 PROCEDURE :: getscale1
96 PROCEDURE :: getscale2
97 PROCEDURE :: setazimuthindex
98 PROCEDURE :: getazimuthindex
99
100! PUBLIC
101 PROCEDURE, PUBLIC :: initgeometry
102 PROCEDURE, PUBLIC :: finalize_base
103 procedure(finalize), DEFERRED :: finalize
104
106 generic, PUBLIC :: scalefactors => scalefactors_0, scalefactors_1, scalefactors_2, &
107 scalefactors_3, scalefactors_4
109 generic, PUBLIC :: radius => radius_0, radius_1, radius_2, radius_3, radius_4
111 generic, PUBLIC :: positionvector => positionvector_0, positionvector_1, positionvector_2, &
112 positionvector_3, positionvector_4
114 generic, PUBLIC :: convert2cartesian => convert2cartesian_coords, convert2cartesian_vectors, &
115 convert2cartesian_coords_1, convert2cartesian_coords_2,&
116 convert2cartesian_coords_3, convert2cartesian_coords_4, &
117 convert2cartesian_vectors_1, convert2cartesian_vectors_2,&
118 convert2cartesian_vectors_3, convert2cartesian_vectors_4
120 generic, PUBLIC :: convert2curvilinear => convert2curvilinear_coords, convert2curvilinear_vectors, &
121 convert2curvilinear_coords_1, convert2curvilinear_coords_2,&
122 convert2curvilinear_coords_3, convert2curvilinear_coords_4,&
123 convert2curvilinear_vectors_1, convert2curvilinear_vectors_2,&
124 convert2curvilinear_vectors_3, convert2curvilinear_vectors_4
125 generic :: setscale => setscale1, setscale2, setscale3
126 generic :: getscale => getscale1, getscale2
127 END TYPE geometry_base
128
129 abstract INTERFACE
130 PURE SUBROUTINE scalefactors_1(this,coords,hx,hy,hz)
131 IMPORT geometry_base
132 IMPLICIT NONE
133 CLASS(geometry_base), INTENT(IN) :: this
134 REAL, DIMENSION(:,:), INTENT(IN) :: coords
135 REAL, DIMENSION(:), INTENT(OUT) :: hx,hy,hz
136 END SUBROUTINE
137 PURE SUBROUTINE scalefactors_2(this,coords,hx,hy,hz)
138 IMPORT geometry_base
139 IMPLICIT NONE
140 CLASS(geometry_base), INTENT(IN) :: this
141 REAL, DIMENSION(:,:,:), INTENT(IN) :: coords
142 REAL, DIMENSION(:,:), INTENT(OUT) :: hx,hy,hz
143 END SUBROUTINE
144 PURE SUBROUTINE scalefactors_3(this,coords,hx,hy,hz)
145 IMPORT geometry_base
146 IMPLICIT NONE
147 CLASS(geometry_base), INTENT(IN) :: this
148 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: coords
149 REAL, DIMENSION(:,:,:), INTENT(OUT) :: hx,hy,hz
150 END SUBROUTINE
151 PURE SUBROUTINE scalefactors_4(this,coords,hx,hy,hz)
152 IMPORT geometry_base
153 IMPLICIT NONE
154 CLASS(geometry_base), INTENT(IN) :: this
155 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: coords
156 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: hx,hy,hz
157 END SUBROUTINE
158 PURE SUBROUTINE radius_1(this,coords,r)
159 IMPORT geometry_base
160 IMPLICIT NONE
161 CLASS(geometry_base), INTENT(IN) :: this
162 REAL, DIMENSION(:,:), INTENT(IN) :: coords
163 REAL, DIMENSION(:), INTENT(OUT) :: r
164 END SUBROUTINE
165 PURE SUBROUTINE radius_2(this,coords,r)
166 IMPORT geometry_base
167 IMPLICIT NONE
168 CLASS(geometry_base), INTENT(IN) :: this
169 REAL, DIMENSION(:,:,:), INTENT(IN) :: coords
170 REAL, DIMENSION(:,:), INTENT(OUT) :: r
171 END SUBROUTINE
172 PURE SUBROUTINE radius_3(this,coords,r)
173 IMPORT geometry_base
174 IMPLICIT NONE
175 CLASS(geometry_base), INTENT(IN) :: this
176 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: coords
177 REAL, DIMENSION(:,:,:), INTENT(OUT) :: r
178 END SUBROUTINE
179 PURE SUBROUTINE radius_4(this,coords,r)
180 IMPORT geometry_base
181 IMPLICIT NONE
182 CLASS(geometry_base), INTENT(IN) :: this
183 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: coords
184 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: r
185 END SUBROUTINE
186 PURE SUBROUTINE positionvector_1(this,coords,posvec)
187 IMPORT geometry_base
188 IMPLICIT NONE
189 CLASS(geometry_base), INTENT(IN) :: this
190 REAL, DIMENSION(:,:), INTENT(IN) :: coords
191 REAL, DIMENSION(:,:), INTENT(OUT) :: posvec
192 END SUBROUTINE
193 PURE SUBROUTINE positionvector_2(this,coords,posvec)
194 IMPORT geometry_base
195 IMPLICIT NONE
196 CLASS(geometry_base), INTENT(IN) :: this
197 REAL, DIMENSION(:,:,:), INTENT(IN) :: coords
198 REAL, DIMENSION(:,:,:), INTENT(OUT) :: posvec
199 END SUBROUTINE
200 PURE SUBROUTINE positionvector_3(this,coords,posvec)
201 IMPORT geometry_base
202 IMPLICIT NONE
203 CLASS(geometry_base), INTENT(IN) :: this
204 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: coords
205 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: posvec
206 END SUBROUTINE
207 PURE SUBROUTINE positionvector_4(this,coords,posvec)
208 IMPORT geometry_base
209 IMPLICIT NONE
210 CLASS(geometry_base), INTENT(IN) :: this
211 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: coords
212 REAL, DIMENSION(:,:,:,:,:), INTENT(OUT) :: posvec
213 END SUBROUTINE
214 PURE SUBROUTINE convert2cartesian_coords_1(this,curv,cart)
215 IMPORT geometry_base
216 IMPLICIT NONE
217 CLASS(geometry_base), INTENT(IN) :: this
218 REAL, DIMENSION(:,:), INTENT(IN) :: curv
219 REAL, DIMENSION(:,:), INTENT(OUT) :: cart
220 END SUBROUTINE
221 PURE SUBROUTINE convert2cartesian_coords_2(this,curv,cart)
222 IMPORT geometry_base
223 IMPLICIT NONE
224 CLASS(geometry_base), INTENT(IN) :: this
225 REAL, DIMENSION(:,:,:), INTENT(IN) :: curv
226 REAL, DIMENSION(:,:,:), INTENT(OUT) :: cart
227 END SUBROUTINE
228 PURE SUBROUTINE convert2cartesian_coords_3(this,curv,cart)
229 IMPORT geometry_base
230 IMPLICIT NONE
231 CLASS(geometry_base), INTENT(IN) :: this
232 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: curv
233 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: cart
234 END SUBROUTINE
235 PURE SUBROUTINE convert2cartesian_coords_4(this,curv,cart)
236 IMPORT geometry_base
237 IMPLICIT NONE
238 CLASS(geometry_base), INTENT(IN) :: this
239 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: curv
240 REAL, DIMENSION(:,:,:,:,:), INTENT(OUT) :: cart
241 END SUBROUTINE
242 PURE SUBROUTINE convert2curvilinear_coords_1(this,cart,curv)
243 IMPORT geometry_base
244 IMPLICIT NONE
245 CLASS(geometry_base), INTENT(IN) :: this
246 REAL, DIMENSION(:,:), INTENT(IN) :: cart
247 REAL, DIMENSION(:,:), INTENT(OUT) :: curv
248 END SUBROUTINE
249 PURE SUBROUTINE convert2curvilinear_coords_2(this,cart,curv)
250 IMPORT geometry_base
251 IMPLICIT NONE
252 CLASS(geometry_base), INTENT(IN) :: this
253 REAL, DIMENSION(:,:,:), INTENT(IN) :: cart
254 REAL, DIMENSION(:,:,:), INTENT(OUT) :: curv
255 END SUBROUTINE
256 PURE SUBROUTINE convert2curvilinear_coords_3(this,cart,curv)
257 IMPORT geometry_base
258 IMPLICIT NONE
259 CLASS(geometry_base), INTENT(IN) :: this
260 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: cart
261 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: curv
262 END SUBROUTINE
263 PURE SUBROUTINE convert2curvilinear_coords_4(this,cart,curv)
264 IMPORT geometry_base
265 IMPLICIT NONE
266 CLASS(geometry_base), INTENT(IN) :: this
267 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: cart
268 REAL, DIMENSION(:,:,:,:,:), INTENT(OUT) :: curv
269 END SUBROUTINE
270 PURE SUBROUTINE convert2cartesian_vectors_1(this,curv,v_curv,v_cart)
271 IMPORT geometry_base
272 IMPLICIT NONE
273 CLASS(geometry_base), INTENT(IN) :: this
274 REAL, DIMENSION(:,:), INTENT(IN) :: curv
275 REAL, DIMENSION(:,:), INTENT(IN) :: v_curv
276 REAL, DIMENSION(:,:), INTENT(OUT) :: v_cart
277 END SUBROUTINE
278 PURE SUBROUTINE convert2cartesian_vectors_2(this,curv,v_curv,v_cart)
279 IMPORT geometry_base
280 IMPLICIT NONE
281 CLASS(geometry_base), INTENT(IN) :: this
282 REAL, DIMENSION(:,:,:), INTENT(IN) :: curv
283 REAL, DIMENSION(:,:,:), INTENT(IN) :: v_curv
284 REAL, DIMENSION(:,:,:), INTENT(OUT) :: v_cart
285 END SUBROUTINE
286 PURE SUBROUTINE convert2cartesian_vectors_3(this,curv,v_curv,v_cart)
287 IMPORT geometry_base
288 IMPLICIT NONE
289 CLASS(geometry_base), INTENT(IN) :: this
290 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: curv
291 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: v_curv
292 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: v_cart
293 END SUBROUTINE
294 PURE SUBROUTINE convert2cartesian_vectors_4(this,curv,v_curv,v_cart)
295 IMPORT geometry_base
296 IMPLICIT NONE
297 CLASS(geometry_base), INTENT(IN) :: this
298 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: curv
299 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: v_curv
300 REAL, DIMENSION(:,:,:,:,:), INTENT(OUT) :: v_cart
301 END SUBROUTINE
302 PURE SUBROUTINE convert2curvilinear_vectors_1(this,curv,v_cart,v_curv)
303 IMPORT geometry_base
304 IMPLICIT NONE
305 CLASS(geometry_base), INTENT(IN) :: this
306 REAL, DIMENSION(:,:), INTENT(IN) :: curv
307 REAL, DIMENSION(:,:), INTENT(IN) :: v_cart
308 REAL, DIMENSION(:,:), INTENT(OUT) :: v_curv
309 END SUBROUTINE
310 PURE SUBROUTINE convert2curvilinear_vectors_2(this,curv,v_cart,v_curv)
311 IMPORT geometry_base
312 IMPLICIT NONE
313 CLASS(geometry_base), INTENT(IN) :: this
314 REAL, DIMENSION(:,:,:), INTENT(IN) :: curv
315 REAL, DIMENSION(:,:,:), INTENT(IN) :: v_cart
316 REAL, DIMENSION(:,:,:), INTENT(OUT) :: v_curv
317 END SUBROUTINE
318 PURE SUBROUTINE convert2curvilinear_vectors_3(this,curv,v_cart,v_curv)
319 IMPORT geometry_base
320 IMPLICIT NONE
321 CLASS(geometry_base), INTENT(IN) :: this
322 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: curv
323 REAL, DIMENSION(:,:,:,:), INTENT(IN) :: v_cart
324 REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: v_curv
325 END SUBROUTINE
326 PURE SUBROUTINE convert2curvilinear_vectors_4(this,curv,v_cart,v_curv)
327 IMPORT geometry_base
328 IMPLICIT NONE
329 CLASS(geometry_base), INTENT(IN) :: this
330 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: curv
331 REAL, DIMENSION(:,:,:,:,:), INTENT(IN) :: v_cart
332 REAL, DIMENSION(:,:,:,:,:), INTENT(OUT) :: v_curv
333 END SUBROUTINE
334 SUBROUTINE finalize(this)
335 IMPORT geometry_base
336 IMPLICIT NONE
337 CLASS(geometry_base), INTENT(INOUT) :: this
338 END SUBROUTINE
339 END INTERFACE
340
344 !--------------------------------------------------------------------------!
345 INTEGER, PARAMETER :: cartesian = 1
346! INTEGER, PARAMETER :: SINHCARTESIAN = 2
347! INTEGER, PARAMETER :: POLAR = 20
348! INTEGER, PARAMETER :: LOGPOLAR = 21
349! INTEGER, PARAMETER :: TANPOLAR = 22
350! INTEGER, PARAMETER :: SINHPOLAR = 23
351! INTEGER, PARAMETER :: SINHTANHPOLAR = 24
352! INTEGER, PARAMETER :: POLYPOLAR = 25
353! INTEGER, PARAMETER :: ELLIPTIC = 27
354! INTEGER, PARAMETER :: BIPOLAR = 29
355 INTEGER, PARAMETER :: cylindrical = 30
356 INTEGER, PARAMETER :: logcylindrical = 31
357 INTEGER, PARAMETER :: tancylindrical = 32
358! INTEGER, PARAMETER :: LNCOSHCYLINDRICAL = 33
359 INTEGER, PARAMETER :: spherical = 40
360 INTEGER, PARAMETER :: logspherical = 41
361 INTEGER, PARAMETER :: spherical_planet = 42
362! INTEGER, PARAMETER :: SINHSPHERICAL = 43
363! INTEGER, PARAMETER :: OBLATE_SPHEROIDAL = 50
364! INTEGER, PARAMETER :: CHANNEL = 60
366 !--------------------------------------------------------------------------!
367 PUBLIC :: &
368 ! types
370 pi, &
371 cartesian, &
374 !--------------------------------------------------------------------------!
375
376CONTAINS
377
379 SUBROUTINE initgeometry(this,gnum,gname,config)
380 IMPLICIT NONE
381 !------------------------------------------------------------------------!
382 CLASS(geometry_base), INTENT(INOUT) :: this
383 TYPE(dict_typ),POINTER :: config
384 INTEGER :: gnum, gt !\todo{same variable}
385 CHARACTER(LEN=*) :: gname
386 REAL :: gs_def,gs_def2,gs_def3
387 !------------------------------------------------------------------------!
388 CHARACTER(LEN=8) :: gs_str
389 !------------------------------------------------------------------------!
390 INTENT(IN) :: gnum,gname
391 !------------------------------------------------------------------------!
392 CALL this%logging_base%InitLogging(gnum,gname)
393
394 CALL getattr(config, "geometry", gt)
395
396 ! check if geometry parameters were given
397 ! and set to defaults if not
398 SELECT CASE(gnum)
400 ! do nothing (no parameters needed)
401 CASE DEFAULT
402 ! geometries with at least one parameter
403 ! gparam defaults to 1 except for CHANNEL
404! IF (gt.EQ.CHANNEL) THEN
405! CALL GetAttr(config, "gparam", gs_def, 0.5)
406! ELSE
407 CALL getattr(config, "gparam", gs_def, 1.0)
408! END IF
409 ! geometries with three parameters
410! SELECT CASE(gt)
411! CASE(SINHPOLAR)
412! ! gp2,gp3 defaults to 0
413! CALL GetAttr(config, "gparam2", gs_def2, 0.0)
414! CALL GetAttr(config, "gparam3", gs_def3, 0.0)
415! CASE(SINHTANHPOLAR,LNCOSHCYLINDRICAL)
416! ! gp2,gp3 defaults to 1
417! CALL GetAttr(config, "gparam2", gs_def2, 1.0)
418! CALL GetAttr(config, "gparam3", gs_def3, 1.0)
419! END SELECT
420 END SELECT
421
422 ! print some information
423 CALL this%Info(" GEOMETRY-> coordinates: " // trim(this%GetName()))
424 SELECT CASE(gt)
426 WRITE (gs_str,'(ES8.1)') this%GetScale()
427 CALL this%Info( " geometry scale: " // trim(gs_str))
428! CASE(LNCOSHCYLINDRICAL,SINHTANHPOLAR)
429! WRITE (gs_str,'(ES8.1)') this%GetScale()
430! CALL this%Info( " geometry scale: " // TRIM(gs_str))
431! WRITE (gs_str,'(ES8.1)') this%GetScale(2)
432! CALL this%Info(" geometry scale: " // TRIM(gs_str))
433! WRITE (gs_str,'(ES8.1)') this%GetScale(3)
434! CALL this%Info(" geometry scale: " // TRIM(gs_str))
435
436 END SELECT
437 END SUBROUTINE initgeometry
438
439 PURE FUNCTION getscale1(this) RESULT(gp)
440 IMPLICIT NONE
441 !------------------------------------------------------------------------!
442 CLASS(geometry_base), INTENT(IN) :: this
443 REAL :: gp
444 !------------------------------------------------------------------------!
445 gp = this%geoparam(1)
446 END FUNCTION getscale1
447
448 PURE FUNCTION getscale2(this,i) RESULT(gp)
449 IMPLICIT NONE
450 !------------------------------------------------------------------------!
451 CLASS(geometry_base), INTENT(IN) :: this
452 INTEGER, INTENT(IN) :: i
453 REAL :: gp
454 !------------------------------------------------------------------------!
455 gp = this%geoparam(i)
456 END FUNCTION getscale2
457
458 PURE SUBROUTINE setscale1(this,gp)
459 IMPLICIT NONE
460 !------------------------------------------------------------------------!
461 CLASS(geometry_base), INTENT(INOUT) :: this
462 REAL, INTENT(IN) :: gp
463 !------------------------------------------------------------------------!
464 this%geoparam(1) = gp
465 END SUBROUTINE setscale1
466
467 PURE SUBROUTINE setscale2(this,gp,gp2)
468 IMPLICIT NONE
469 !------------------------------------------------------------------------!
470 CLASS(geometry_base), INTENT(INOUT) :: this
471 REAL, INTENT(IN) :: gp,gp2
472 !------------------------------------------------------------------------!
473 this%geoparam(1) = gp
474 this%geoparam(2) = gp2
475 END SUBROUTINE setscale2
476
477 PURE SUBROUTINE setscale3(this,gp,gp2,gp3)
478 IMPLICIT NONE
479 !------------------------------------------------------------------------!
480 CLASS(geometry_base), INTENT(INOUT) :: this
481 REAL, INTENT(IN) :: gp,gp2,gp3
482 !------------------------------------------------------------------------!
483 this%geoparam(1) = gp
484 this%geoparam(2) = gp2
485 this%geoparam(3) = gp3
486 END SUBROUTINE setscale3
487
488
499 PURE SUBROUTINE scalefactors_0(this,coords,hx,hy,hz)
500 IMPLICIT NONE
501 !------------------------------------------------------------------------!
502 CLASS(geometry_base), INTENT(IN) :: this
503 TYPE(marray_cellvector), INTENT(IN) :: coords
504 TYPE(marray_cellscalar), INTENT(INOUT) :: hx,hy,hz
505 !------------------------------------------------------------------------!
506 CALL this%ScaleFactors(coords%data2d(:,:),hx%data1d(:),hy%data1d(:),hz%data1d(:))
507 END SUBROUTINE scalefactors_0
508
515 PURE SUBROUTINE radius_0(this,coords,radius)
516 IMPLICIT NONE
517 !------------------------------------------------------------------------!
518 CLASS(geometry_base), INTENT(IN) :: this
519 TYPE(marray_cellvector), INTENT(IN) :: coords
520 TYPE(marray_cellscalar), INTENT(INOUT) :: radius
521 !------------------------------------------------------------------------!
522 CALL this%Radius(coords%data2d(:,:),radius%data1d(:))
523 END SUBROUTINE radius_0
524
531 PURE SUBROUTINE positionvector_0(this,coords,posvec)
532 IMPLICIT NONE
533 !------------------------------------------------------------------------!
534 CLASS(geometry_base), INTENT(IN) :: this
535 TYPE(marray_cellvector), INTENT(IN) :: coords
536 TYPE(marray_cellvector), INTENT(INOUT) :: posvec
537 !------------------------------------------------------------------------!
538 CALL this%PositionVector(coords%data2d(:,:),posvec%data2d(:,:))
539 END SUBROUTINE positionvector_0
540
545 PURE SUBROUTINE convert2cartesian_coords(this,curv,cart)
546 IMPLICIT NONE
547 !------------------------------------------------------------------------!
548 CLASS(geometry_base), INTENT(IN) :: this
549 TYPE(marray_cellvector), INTENT(IN) :: curv
550 TYPE(marray_cellvector), INTENT(INOUT) :: cart
551 !------------------------------------------------------------------------!
552 CALL this%Convert2Cartesian(curv%data2d(:,:),cart%data2d(:,:))
553 END SUBROUTINE convert2cartesian_coords
554
559 PURE SUBROUTINE convert2curvilinear_coords(this,cart,curv)
560 IMPLICIT NONE
561 !------------------------------------------------------------------------!
562 CLASS(geometry_base), INTENT(IN) :: this
563 TYPE(marray_cellvector), INTENT(IN) :: cart
564 TYPE(marray_cellvector), INTENT(INOUT) :: curv
565 !------------------------------------------------------------------------!
566 CALL this%Convert2Curvilinear(cart%data2d(:,:),curv%data2d(:,:))
567 END SUBROUTINE convert2curvilinear_coords
568
573 PURE SUBROUTINE convert2cartesian_vectors(this,curv,v_curv,v_cart)
574 IMPLICIT NONE
575 !------------------------------------------------------------------------!
576 CLASS(geometry_base), INTENT(IN) :: this
577 TYPE(marray_cellvector), INTENT(IN) :: curv
578 TYPE(marray_cellvector), INTENT(IN) :: v_curv
579 TYPE(marray_cellvector), INTENT(INOUT) :: v_cart
580 !------------------------------------------------------------------------!
581 CALL this%Convert2Cartesian(curv%data2d(:,:),v_curv%data2d(:,:),v_cart%data2d(:,:))
582 END SUBROUTINE convert2cartesian_vectors
583
584
586 PURE SUBROUTINE convert2curvilinear_vectors(this,curv,v_cart,v_curv)
587 IMPLICIT NONE
588 !------------------------------------------------------------------------!
589 CLASS(geometry_base), INTENT(IN) :: this
590 TYPE(marray_cellvector), INTENT(IN) :: curv
591 TYPE(marray_cellvector), INTENT(IN) :: v_cart
592 TYPE(marray_cellvector), INTENT(INOUT) :: v_curv
593 !------------------------------------------------------------------------!
594 CALL this%Convert2Curvilinear(curv%data2d(:,:),v_cart%data2d(:,:),v_curv%data2d(:,:))
595 END SUBROUTINE convert2curvilinear_vectors
596
600 PURE SUBROUTINE setazimuthindex(this,idx)
601 IMPLICIT NONE
602 !------------------------------------------------------------------------!
603 CLASS(geometry_base), INTENT(INOUT) :: this
604 INTEGER, INTENT(IN) :: idx
605 !------------------------------------------------------------------------!
606 this%azimuthIndex = idx
607 END SUBROUTINE setazimuthindex
608
612 PURE FUNCTION getazimuthindex(this) RESULT(idx)
613 IMPLICIT NONE
614 !------------------------------------------------------------------------!
615 CLASS(geometry_base), INTENT(IN) :: this
616 INTEGER :: idx
617 !------------------------------------------------------------------------!
618 idx = this%azimuthIndex
619 END FUNCTION getazimuthindex
620
622 SUBROUTINE finalize_base(this)
623 IMPLICIT NONE
624 !------------------------------------------------------------------------!
625 CLASS(geometry_base), INTENT(INOUT) :: this
626 !------------------------------------------------------------------------!
627 IF (.NOT.this%Initialized()) &
628 CALL this%Error("CloseGeometry","not initialized")
629 !CALL this%logging_base%Finalize()
630 END SUBROUTINE finalize_base
631
632END MODULE geometry_base_mod
Dictionary for generic data types.
Definition: common_dict.f90:61
type(logging_base), save this
base class for geometrical properties
pure integer function getazimuthindex(this)
returns the coordinate index of the azimuthal angle
real, parameter, public pi
integer, parameter, public spherical
pure subroutine convert2cartesian_vectors(this, curv, v_curv, v_cart)
Convert curvilinear vector components to cartesian vector components.
integer, parameter, public tancylindrical
pure subroutine positionvector_0(this, coords, posvec)
compute position vector components for all cell positions
integer, parameter, public cartesian
integer, parameter, public cylindrical
subroutine finalize_base(this)
Destructor of generic geometry module.
integer, parameter, public logcylindrical
subroutine initgeometry(this, gnum, gname, config)
Constructor of generic geometry module.
pure subroutine setscale2(this, gp, gp2)
pure subroutine scalefactors_0(this, coords, hx, hy, hz)
Compute scale factors.
pure subroutine radius_0(this, coords, radius)
Compute radial distances to the origin.
pure subroutine convert2curvilinear_coords(this, cart, curv)
Convert cartesian to curvilinear coordinates.
pure real function getscale2(this, i)
integer, parameter, public spherical_planet
pure subroutine convert2curvilinear_vectors(this, curv, v_cart, v_curv)
Convert cartesian vector components to curvilinear vector components.
pure subroutine convert2cartesian_coords(this, curv, cart)
Convert curvilinear to cartesian coordinates.
integer, parameter, public logspherical
pure real function getscale1(this)
pure subroutine setazimuthindex(this, idx)
sets the coordinate index of the azimuthal angle
pure subroutine setscale1(this, gp)
pure subroutine setscale3(this, gp, gp2, gp3)
Basic fosite module.
subroutine finalize(this)
Destructor of logging_base class.
derived mesh array class for scalar cell data
derived mesh array class for vector cell data
common data structure