geomechanics_regression.F90 coverage: 85.71 %func 47.27 %block
1) module Geomechanics_Regression_module
2)
3) use Output_Aux_module
4)
5) use PFLOTRAN_Constants_module
6)
7) implicit none
8)
9) private
10)
11) #include "petsc/finclude/petscsys.h"
12) #include "petsc/finclude/petscvec.h"
13) #include "petsc/finclude/petscvec.h90"
14)
15) type, public :: geomechanics_regression_type
16) type(geomechanics_regression_variable_type), pointer :: variable_list
17) PetscInt, pointer :: natural_vertex_ids(:)
18) PetscInt :: num_vertices_per_process
19) PetscInt, pointer :: vertices_per_process_natural_ids(:)
20) Vec :: natural_vertex_id_vec
21) Vec :: vertices_per_process_vec
22) VecScatter :: scatter_natural_vertex_id_gtos
23) VecScatter :: scatter_vertices_per_process_gtos
24) type(geomechanics_regression_type), pointer :: next
25) end type geomechanics_regression_type
26)
27) type, public :: geomechanics_regression_variable_type
28) character(len=MAXSTRINGLENGTH) :: name
29) type(geomechanics_regression_variable_type), pointer :: next
30) end type geomechanics_regression_variable_type
31)
32) public :: GeomechanicsRegressionRead, &
33) GeomechanicsRegressionCreateMapping, &
34) GeomechanicsRegressionOutput, &
35) GeomechanicsRegressionDestroy
36)
37) contains
38)
39) ! ************************************************************************** !
40)
41) function GeomechanicsRegressionCreate()
42) !
43) ! Creates a geomechanics regression object
44) !
45) ! Author: Satish Karra
46) ! Date: 06/01/2016
47) !
48)
49) implicit none
50)
51) type(geomechanics_regression_type), pointer :: GeomechanicsRegressionCreate
52)
53) type(geomechanics_regression_type), pointer :: geomechanics_regression
54)
55) allocate(geomechanics_regression)
56) nullify(geomechanics_regression%variable_list)
57) nullify(geomechanics_regression%natural_vertex_ids)
58) geomechanics_regression%num_vertices_per_process = 0
59) nullify(geomechanics_regression%vertices_per_process_natural_ids)
60) geomechanics_regression%natural_vertex_id_vec = 0
61) geomechanics_regression%vertices_per_process_vec = 0
62) geomechanics_regression%scatter_natural_vertex_id_gtos = 0
63) geomechanics_regression%scatter_vertices_per_process_gtos = 0
64) nullify(geomechanics_regression%next)
65) GeomechanicsRegressionCreate => geomechanics_regression
66)
67) end function GeomechanicsRegressionCreate
68)
69) ! ************************************************************************** !
70)
71) function GeomechanicsRegressionVariableCreate()
72) !
73) ! Creates a geomechanics_regression variable object
74) !
75) ! Author: Satish Karra
76) ! Date: 06/01/2016
77) !
78)
79) implicit none
80)
81) type(geomechanics_regression_variable_type), pointer :: GeomechanicsRegressionVariableCreate
82)
83) type(geomechanics_regression_variable_type), pointer :: geomechanics_regression_variable
84)
85) allocate(geomechanics_regression_variable)
86) geomechanics_regression_variable%name = ''
87) nullify(geomechanics_regression_variable%next)
88) GeomechanicsRegressionVariableCreate => geomechanics_regression_variable
89)
90) end function GeomechanicsRegressionVariableCreate
91)
92) ! ************************************************************************** !
93)
94) subroutine GeomechanicsRegressionRead(geomechanics_regression,input,option)
95) !
96) ! Reads in contents of a geomechanics_regression card
97) !
98) ! Author: Satish Karra
99) ! Date: 06/22/2016
100) !
101)
102) use Option_module
103) use Input_Aux_module
104) use String_module
105) use Utility_module
106)
107) implicit none
108)
109) type(geomechanics_regression_type), pointer :: geomechanics_regression
110) type(input_type), pointer :: input
111) type(option_type) :: option
112)
113) character(len=MAXWORDLENGTH) :: keyword, word
114) type(geomechanics_regression_variable_type), pointer :: cur_variable, new_variable
115) PetscInt :: count, max_vertices
116) PetscInt, pointer :: int_array(:)
117) PetscErrorCode :: ierr
118)
119) geomechanics_regression => GeomechanicsRegressionCreate()
120)
121) input%ierr = 0
122) do
123)
124) call InputReadPflotranString(input,option)
125)
126) if (InputCheckExit(input,option)) exit
127)
128) call InputReadWord(input,option,keyword,PETSC_TRUE)
129) call InputErrorMsg(input,option,'keyword','GEOMECHANICS_REGRESSION')
130) call StringToUpper(keyword)
131)
132) select case(trim(keyword))
133)
134) case('VARIABLES')
135) count = 0
136) do
137) call InputReadPflotranString(input,option)
138) if (InputCheckExit(input,option)) exit
139)
140) call InputReadWord(input,option,word,PETSC_TRUE)
141) call InputErrorMsg(input,option,'variable','GEOMECHANICS_REGRESSION,VARIABLES')
142) call StringToUpper(word)
143) new_variable => GeomechanicsRegressionVariableCreate()
144) new_variable%name = word
145) if (.not.associated(geomechanics_regression%variable_list)) then
146) geomechanics_regression%variable_list => new_variable
147) else
148) cur_variable%next => new_variable
149) endif
150) cur_variable => new_variable
151) enddo
152) case('VERTICES')
153) max_vertices = 100
154) allocate(int_array(max_vertices))
155) count = 0
156) do
157) call InputReadPflotranString(input,option)
158) if (InputCheckExit(input,option)) exit
159)
160) count = count + 1
161) if (count > max_vertices) then
162) call reallocateIntArray(int_array,max_vertices)
163) endif
164) call InputReadInt(input,option,int_array(count))
165) call InputErrorMsg(input,option,'natural vertex id','GEOMECHANICS_REGRESSION,VERTICES')
166) enddo
167) allocate(geomechanics_regression%natural_vertex_ids(count))
168) geomechanics_regression%natural_vertex_ids = int_array(1:count)
169) call PetscSortInt(count,geomechanics_regression%natural_vertex_ids, &
170) ierr);CHKERRQ(ierr)
171) deallocate(int_array)
172) case('VERTICES_PER_PROCESS')
173) call InputReadInt(input,option,geomechanics_regression%num_vertices_per_process)
174) call InputErrorMsg(input,option,'num vertices per process','GEOMECHANICS_REGRESSION')
175) case default
176) call InputKeywordUnrecognized(keyword,'GEOMECHANICS_REGRESSION',option)
177) end select
178)
179) enddo
180)
181) end subroutine GeomechanicsRegressionRead
182)
183) ! ************************************************************************** !
184)
185) subroutine GeomechanicsRegressionCreateMapping(geomechanics_regression, &
186) geomechanics_realization)
187) !
188) ! Creates mapping between a natural mpi vec and a
189) ! sequential vec on io_rank
190) !
191) ! Author: Satish Karra
192) ! Date: 06/22/2016
193) !
194)
195) use Option_module
196) use Geomechanics_Realization_class
197) use Geomechanics_Grid_Aux_module
198) use Geomechanics_Discretization_module
199)
200) implicit none
201)
202) #include "petsc/finclude/petscis.h"
203) #include "petsc/finclude/petscis.h90"
204) #include "petsc/finclude/petscviewer.h"
205)
206) type(geomechanics_regression_type), pointer :: geomechanics_regression
207) class(realization_geomech_type) :: geomechanics_realization
208)
209) IS :: is_petsc
210) PetscInt, allocatable :: int_array(:)
211) PetscInt :: i, upper_bound, lower_bound, count, temp_int
212) PetscInt :: local_id
213) PetscReal, pointer :: vec_ptr(:)
214) Vec :: temp_vec
215) VecScatter :: temp_scatter
216) IS :: temp_is
217) PetscViewer :: viewer
218) PetscErrorCode :: ierr
219)
220) type(geomech_grid_type), pointer :: grid
221) type(option_type), pointer :: option
222)
223) if (.not.associated(geomechanics_regression)) return
224)
225) grid => geomechanics_realization%geomech_patch%geomech_grid
226) option => geomechanics_realization%option
227)
228) ! natural vertex ids
229) if (associated(geomechanics_regression%natural_vertex_ids)) then
230) ! ensure that natural ids are within problem domain
231) if (maxval(geomechanics_regression%natural_vertex_ids) > &
232) grid%nmax_node) then
233) option%io_buffer = 'Natural IDs outside geomechanics domain ' // &
234) 'requested for geomechanics regression output. ' // &
235) 'Removing non-existent IDs.'
236) call printWrnMsg(option)
237) count = 0
238) allocate(int_array(size(geomechanics_regression%natural_vertex_ids)))
239) int_array = 0
240) do i = 1, size(geomechanics_regression%natural_vertex_ids)
241) if (geomechanics_regression%natural_vertex_ids(i) &
242) <= grid%nmax_node) then
243) count = count + 1
244) int_array(count) = geomechanics_regression%natural_vertex_ids(i)
245) endif
246) enddo
247) ! reallocate array
248) deallocate(geomechanics_regression%natural_vertex_ids)
249) allocate(geomechanics_regression%natural_vertex_ids(count))
250) !geh: Since natural_vertex_ids and int_array may now be of different sizes,
251) ! we need to be explicit about the values to copy. gfortran has
252) ! issues with this while Intel figures it out. Better to be explicit.
253) geomechanics_regression%natural_vertex_ids = int_array(1:count)
254) deallocate(int_array)
255) endif
256) ! Create a local vector on io_rank and scatter the nodal data
257) ! to that local vector
258) call VecCreate(PETSC_COMM_SELF, &
259) geomechanics_regression%natural_vertex_id_vec, &
260) ierr);CHKERRQ(ierr)
261) if (option%myrank == option%io_rank) then
262) call VecSetSizes(geomechanics_regression%natural_vertex_id_vec, &
263) size(geomechanics_regression%natural_vertex_ids), &
264) PETSC_DECIDE,ierr);CHKERRQ(ierr)
265) else
266) call VecSetSizes(geomechanics_regression%natural_vertex_id_vec, &
267) ZERO_INTEGER, &
268) PETSC_DECIDE,ierr);CHKERRQ(ierr)
269) endif
270) call VecSetFromOptions(geomechanics_regression%natural_vertex_id_vec, &
271) ierr);CHKERRQ(ierr)
272)
273) if (option%myrank == option%io_rank) then
274) count = size(geomechanics_regression%natural_vertex_ids)
275) ! determine how many of the natural vertex ids are local
276) allocate(int_array(count))
277) int_array = geomechanics_regression%natural_vertex_ids
278) ! convert to zero based
279) int_array = int_array - 1
280) else
281) count = 0
282) allocate(int_array(count))
283) endif
284) call GeomechDiscretAOApplicationToPetsc(geomechanics_realization% &
285) geomech_discretization, &
286) int_array)
287)
288) ! create IS for global petsc vertex ids
289) call ISCreateGeneral(option%mycomm,count,int_array,PETSC_COPY_VALUES, &
290) is_petsc,ierr);CHKERRQ(ierr)
291) deallocate(int_array)
292)
293) #ifdef GEOMECHANICS_REGRESSION_DEBUG
294) call PetscViewerASCIIOpen(option%mycomm, &
295) 'geomech_is_petsc_natural_vertex_id.out', &
296) viewer,ierr);CHKERRQ(ierr)
297) call ISView(is_petsc,viewer,ierr);CHKERRQ(ierr)
298) call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
299) #endif
300)
301) ! create scatter context
302) call VecScatterCreate(geomechanics_realization%geomech_field% &
303) press,is_petsc, &
304) geomechanics_regression%natural_vertex_id_vec, &
305) PETSC_NULL_OBJECT, &
306) geomechanics_regression% &
307) scatter_natural_vertex_id_gtos, &
308) ierr);CHKERRQ(ierr)
309)
310) call ISDestroy(is_petsc,ierr);CHKERRQ(ierr)
311)
312) #ifdef GEOMECHANICS_REGRESSION_DEBUG
313) call PetscViewerASCIIOpen(option%mycomm, &
314) 'geomechanics_regression_scatter_nat_vertex_ids.out', &
315) viewer, &
316) ierr);CHKERRQ(ierr)
317) call VecScatterView(geomechanics_regression% &
318) scatter_natural_vertex_id_gtos,viewer, &
319) ierr);CHKERRQ(ierr)
320) call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
321) #endif
322)
323) endif
324)
325) if (geomechanics_regression%num_vertices_per_process > 0) then
326) ! determine minimum number of vertices per process
327) i = grid%nlmax_node
328) call MPI_Allreduce(i,count,ONE_INTEGER_MPI,MPIU_INTEGER,MPI_MIN, &
329) option%mycomm,ierr)
330) if (count < geomechanics_regression%num_vertices_per_process) then
331) option%io_buffer = 'Number of vertices per process for ' // &
332) 'GeomechanicsRegression exceeds minimum number of vertices ' // &
333) 'per process. Truncating.'
334) call printMsg(option)
335) geomechanics_regression%num_vertices_per_process = count
336) endif
337)
338) ! vertices ids per processor
339) call VecCreate(PETSC_COMM_SELF, &
340) geomechanics_regression%vertices_per_process_vec, &
341) ierr);CHKERRQ(ierr)
342) if (option%myrank == option%io_rank) then
343) call VecSetSizes(geomechanics_regression%vertices_per_process_vec, &
344) geomechanics_regression%num_vertices_per_process* &
345) option%mycommsize, &
346) PETSC_DECIDE,ierr);CHKERRQ(ierr)
347) else
348) call VecSetSizes(geomechanics_regression% &
349) vertices_per_process_vec,ZERO_INTEGER, &
350) PETSC_DECIDE,ierr);CHKERRQ(ierr)
351) endif
352) call VecSetFromOptions(geomechanics_regression%vertices_per_process_vec, &
353) ierr);CHKERRQ(ierr)
354)
355) ! create temporary vec to transfer down ids of vertices
356) call VecCreate(option%mycomm,temp_vec,ierr);CHKERRQ(ierr)
357) call VecSetSizes(temp_vec, &
358) geomechanics_regression%num_vertices_per_process, &
359) PETSC_DECIDE,ierr);CHKERRQ(ierr)
360) call VecSetFromOptions(temp_vec,ierr);CHKERRQ(ierr)
361)
362) ! calculate interval
363) call VecGetArrayF90(temp_vec,vec_ptr,ierr);CHKERRQ(ierr)
364) temp_int = grid%nlmax_node / &
365) geomechanics_regression%num_vertices_per_process
366) do i = 1, geomechanics_regression%num_vertices_per_process
367) vec_ptr(i) = temp_int*(i-1) + 1 + grid%global_offset
368) enddo
369) call VecRestoreArrayF90(temp_vec,vec_ptr,ierr);CHKERRQ(ierr)
370)
371) ! create temporary scatter to transfer values to io_rank
372) if (option%myrank == option%io_rank) then
373) count = option%mycommsize* &
374) geomechanics_regression%num_vertices_per_process
375) ! determine how many of the natural vertex ids are local
376) allocate(int_array(count))
377) do i = 1, count
378) int_array(i) = i
379) enddo
380) ! convert to zero based
381) int_array = int_array - 1
382) else
383) count = 0
384) allocate(int_array(count))
385) endif
386) call ISCreateGeneral(option%mycomm,count, &
387) int_array,PETSC_COPY_VALUES,temp_is, &
388) ierr);CHKERRQ(ierr)
389)
390) call VecScatterCreate(temp_vec,temp_is, &
391) geomechanics_regression%vertices_per_process_vec, &
392) PETSC_NULL_OBJECT, &
393) temp_scatter,ierr);CHKERRQ(ierr)
394) call ISDestroy(temp_is,ierr);CHKERRQ(ierr)
395)
396) ! scatter ids to io_rank
397) call VecScatterBegin(temp_scatter,temp_vec, &
398) geomechanics_regression%vertices_per_process_vec, &
399) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
400) call VecScatterEnd(temp_scatter,temp_vec, &
401) geomechanics_regression%vertices_per_process_vec, &
402) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
403) call VecScatterDestroy(temp_scatter,ierr);CHKERRQ(ierr)
404) call VecDestroy(temp_vec,ierr);CHKERRQ(ierr)
405)
406) ! transfer vertex ids into array for creating new scatter
407) if (option%myrank == option%io_rank) then
408) count = option%mycommsize* &
409) geomechanics_regression%num_vertices_per_process
410) call VecGetArrayF90(geomechanics_regression%vertices_per_process_vec, &
411) vec_ptr,ierr);CHKERRQ(ierr)
412) do i = 1, count
413) int_array(i) = int(vec_ptr(i)+0.1d0) ! tolerance to ensure int value
414) enddo
415) call VecRestoreArrayF90(geomechanics_regression% &
416) vertices_per_process_vec,vec_ptr, &
417) ierr);CHKERRQ(ierr)
418) ! convert to zero based
419) int_array = int_array - 1
420) endif
421)
422) call ISCreateGeneral(option%mycomm,count, &
423) int_array,PETSC_COPY_VALUES,is_petsc, &
424) ierr);CHKERRQ(ierr)
425) deallocate(int_array)
426)
427) #ifdef GEOMECHANICS_REGRESSION_DEBUG
428) call PetscViewerASCIIOpen(option%mycomm, &
429) 'geomech_is_petsc_vertices_per_process.out', &
430) viewer,ierr);CHKERRQ(ierr)
431) call ISView(is_petsc,viewer,ierr);CHKERRQ(ierr)
432) call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
433) #endif
434)
435) call VecScatterCreate(geomechanics_realization%geomech_field% &
436) press,is_petsc, &
437) geomechanics_regression%vertices_per_process_vec, &
438) PETSC_NULL_OBJECT, &
439) geomechanics_regression% &
440) scatter_vertices_per_process_gtos, &
441) ierr);CHKERRQ(ierr)
442) call ISDestroy(is_petsc,ierr);CHKERRQ(ierr)
443)
444) #ifdef GEOMECHANICS_REGRESSION_DEBUG
445) call PetscViewerASCIIOpen(option%mycomm, &
446) 'geomechanics_regression_scatter_vertices_per_process.out', &
447) viewer,ierr);CHKERRQ(ierr)
448) call VecScatterView(geomechanics_regression% &
449) scatter_vertices_per_process_gtos,viewer, &
450) ierr);CHKERRQ(ierr)
451) call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
452) #endif
453)
454) ! fill in natural ids of these vertices on the io_rank
455) if (option%myrank == option%io_rank) then
456) allocate(geomechanics_regression%vertices_per_process_natural_ids( &
457) geomechanics_regression%num_vertices_per_process* &
458) option%mycommsize))
459) endif
460)
461) call VecGetArrayF90(geomechanics_realization%geomech_field%press,vec_ptr, &
462) ierr);CHKERRQ(ierr)
463) do local_id = 1, grid%nlmax_node
464) vec_ptr(local_id) = grid%nG2A(grid%nL2G(local_id))
465) enddo
466) call VecRestoreArrayF90(geomechanics_realization%geomech_field%press, &
467) vec_ptr,ierr);CHKERRQ(ierr)
468)
469) call VecScatterBegin(geomechanics_regression% &
470) scatter_vertices_per_process_gtos, &
471) geomechanics_realization%geomech_field%press, &
472) geomechanics_regression%vertices_per_process_vec, &
473) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
474) call VecScatterEnd(geomechanics_regression% &
475) scatter_vertices_per_process_gtos, &
476) geomechanics_realization%geomech_field%press, &
477) geomechanics_regression%vertices_per_process_vec, &
478) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
479)
480) if (option%myrank == option%io_rank) then
481) call VecGetArrayF90(geomechanics_regression% &
482) vertices_per_process_vec,vec_ptr, &
483) ierr);CHKERRQ(ierr)
484) geomechanics_regression%vertices_per_process_natural_ids(:) = &
485) int(vec_ptr(:)+0.1)
486) call VecRestoreArrayF90(geomechanics_regression% &
487) vertices_per_process_vec,vec_ptr, &
488) ierr);CHKERRQ(ierr)
489) endif
490)
491) endif
492)
493) end subroutine GeomechanicsRegressionCreateMapping
494)
495) ! ************************************************************************** !
496)
497) subroutine GeomechanicsRegressionOutput(geomechanics_regression, &
498) geomechanics_realization, &
499) geomechanics_timestepper)
500) !
501) ! Prints geomechanics regression output through the io_rank
502) !
503) ! Author: Satish Karra
504) ! Date: 06/22/2016
505) !
506)
507) use Geomechanics_Realization_class
508) use Timestepper_Steady_class
509) use Option_module
510) use Geomechanics_Discretization_module
511) use Output_Geomechanics_module, only : OutputGeomechGetVarFromArray
512)
513) implicit none
514)
515) type(geomechanics_regression_type), pointer :: geomechanics_regression
516) class(realization_geomech_type) :: geomechanics_realization
517) class(timestepper_steady_type), pointer :: geomechanics_timestepper
518) ! these must be pointers as they can be null
519) character(len=MAXSTRINGLENGTH) :: string
520) Vec :: global_vec
521) PetscInt :: ivar, isubvar
522) type(option_type), pointer :: option
523) type(output_variable_type), pointer :: cur_variable
524) PetscReal, pointer :: vec_ptr(:), y_ptr(:), z_ptr(:)
525) PetscInt :: i
526) PetscInt :: iphase
527) PetscReal :: r_norm, x_norm
528) PetscReal :: max, min, mean
529) PetscErrorCode :: ierr
530)
531) if (.not.associated(geomechanics_regression)) return
532)
533) option => geomechanics_realization%option
534)
535) if (option%myrank == option%io_rank) then
536) string = trim(option%global_prefix) // &
537) trim(option%group_prefix) // &
538) '.regression'
539) option%io_buffer = '--> write geomechanics_regression output file: ' // trim(string)
540) call printMsg(option)
541) open(unit=OUTPUT_UNIT,file=string,action="write")
542) endif
543)
544) call GeomechDiscretizationCreateVector(geomechanics_realization% &
545) geomech_discretization,ONEDOF, &
546) global_vec,GLOBAL,option)
547) cur_variable => geomechanics_realization%output_option% &
548) output_snap_variable_list%first
549) do
550) if (.not.associated(cur_variable)) exit
551)
552) ivar = cur_variable%ivar
553) isubvar = cur_variable%isubvar
554)
555) call OutputGeomechGetVarFromArray(geomechanics_realization,global_vec,ivar,isubvar)
556) call VecMax(global_vec,PETSC_NULL_INTEGER,max,ierr);CHKERRQ(ierr)
557) call VecMin(global_vec,PETSC_NULL_INTEGER,min,ierr);CHKERRQ(ierr)
558) call VecSum(global_vec,mean,ierr);CHKERRQ(ierr)
559) mean = mean / geomechanics_realization%geomech_patch%geomech_grid%nmax_node
560)
561) ! list of natural ids
562) if (associated(geomechanics_regression%natural_vertex_ids)) then
563) call VecScatterBegin(geomechanics_regression%scatter_natural_vertex_id_gtos, &
564) global_vec, &
565) geomechanics_regression%natural_vertex_id_vec, &
566) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
567) call VecScatterEnd(geomechanics_regression%scatter_natural_vertex_id_gtos, &
568) global_vec, &
569) geomechanics_regression%natural_vertex_id_vec, &
570) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
571) endif
572) if (geomechanics_regression%num_vertices_per_process > 0) then
573) ! vertices per process
574) call VecScatterBegin(geomechanics_regression%scatter_vertices_per_process_gtos, &
575) global_vec, &
576) geomechanics_regression%vertices_per_process_vec, &
577) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
578) call VecScatterEnd(geomechanics_regression%scatter_vertices_per_process_gtos, &
579) global_vec, &
580) geomechanics_regression%vertices_per_process_vec, &
581) INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
582) endif
583)
584) 100 format(i9,': ',es21.13)
585) 101 format(i9,': ',i9)
586)
587) if (option%myrank == option%io_rank) then
588) string = OutputVariableToCategoryString(cur_variable%icategory)
589) write(OUTPUT_UNIT,'(''-- '',a,'': '',a,'' --'')') &
590) trim(string), trim(cur_variable%name)
591)
592) ! max, min, mean
593) if (cur_variable%iformat == 0) then
594) write(OUTPUT_UNIT,'(6x,''Max: '',es21.13)') max
595) write(OUTPUT_UNIT,'(6x,''Min: '',es21.13)') min
596) else
597) write(OUTPUT_UNIT,'(6x,''Max: '',i9)') int(max)
598) write(OUTPUT_UNIT,'(6x,''Min: '',i9)') int(min)
599) endif
600) write(OUTPUT_UNIT,'(5x,''Mean: '',es21.13)') mean
601)
602) ! natural vertex ids
603) if (associated(geomechanics_regression%natural_vertex_ids)) then
604) if (size(geomechanics_regression%natural_vertex_ids) > 0) then
605) call VecGetArrayF90(geomechanics_regression%natural_vertex_id_vec,vec_ptr, &
606) ierr);CHKERRQ(ierr)
607) if (cur_variable%iformat == 0) then
608) do i = 1, size(geomechanics_regression%natural_vertex_ids)
609) write(OUTPUT_UNIT,100) &
610) geomechanics_regression%natural_vertex_ids(i),vec_ptr(i)
611) enddo
612) else
613) do i = 1, size(geomechanics_regression%natural_vertex_ids)
614) write(OUTPUT_UNIT,101) &
615) geomechanics_regression%natural_vertex_ids(i),int(vec_ptr(i))
616) enddo
617) endif
618) call VecRestoreArrayF90(geomechanics_regression%natural_vertex_id_vec,vec_ptr, &
619) ierr);CHKERRQ(ierr)
620) endif
621) endif
622)
623) ! vertex ids per process
624) if (geomechanics_regression%num_vertices_per_process > 0) then
625) call VecGetArrayF90(geomechanics_regression%vertices_per_process_vec,vec_ptr, &
626) ierr);CHKERRQ(ierr)
627) if (cur_variable%iformat == 0) then
628) do i = 1, geomechanics_regression%num_vertices_per_process*option%mycommsize
629) write(OUTPUT_UNIT,100) &
630) geomechanics_regression%vertices_per_process_natural_ids(i),vec_ptr(i)
631) enddo
632) else
633) do i = 1, geomechanics_regression%num_vertices_per_process*option%mycommsize
634) write(OUTPUT_UNIT,101) &
635) geomechanics_regression%vertices_per_process_natural_ids(i),int(vec_ptr(i))
636) enddo
637) endif
638) call VecRestoreArrayF90(geomechanics_regression%vertices_per_process_vec,vec_ptr, &
639) ierr);CHKERRQ(ierr)
640) endif
641) endif
642)
643) cur_variable => cur_variable%next
644) enddo
645)
646) 102 format(i12)
647) 103 format(es21.13)
648)
649) ! timestep, newton iteration, solver iteration output
650) if (associated(geomechanics_timestepper)) then
651) call VecNorm(geomechanics_realization%geomech_field%disp_xx, &
652) NORM_2,x_norm,ierr);CHKERRQ(ierr)
653) call VecNorm(geomechanics_realization%geomech_field%disp_r, &
654) NORM_2,r_norm,ierr);CHKERRQ(ierr)
655) if (option%myrank == option%io_rank) then
656) write(OUTPUT_UNIT,'(''-- SOLUTION: Geomechanics --'')')
657) write(OUTPUT_UNIT,'('' Time (seconds): '',es21.13)') &
658) geomechanics_timestepper%cumulative_solver_time
659) write(OUTPUT_UNIT,'('' Time Steps: '',i12)') geomechanics_timestepper%steps
660) write(OUTPUT_UNIT,'('' Newton Iterations: '',i12)') &
661) geomechanics_timestepper%cumulative_newton_iterations
662) write(OUTPUT_UNIT,'('' Solver Iterations: '',i12)') &
663) geomechanics_timestepper%cumulative_linear_iterations
664) write(OUTPUT_UNIT,'('' Time Step Cuts: '',i12)') &
665) geomechanics_timestepper%cumulative_time_step_cuts
666) write(OUTPUT_UNIT,'('' Solution 2-Norm: '',es21.13)') x_norm
667) write(OUTPUT_UNIT,'('' Residual 2-Norm: '',es21.13)') r_norm
668) endif
669) endif
670)
671) close(OUTPUT_UNIT)
672)
673) end subroutine GeomechanicsRegressionOutput
674)
675) ! ************************************************************************** !
676)
677) recursive subroutine GeomechanicsRegressionVariableDestroy( &
678) geomechanics_regression_variable)
679) !
680) ! Destroys a geomechanics regression variable object
681) !
682) ! Author: Satish Karra
683) ! Date: 06/22/2016
684) !
685)
686) implicit none
687)
688) type(geomechanics_regression_variable_type), pointer :: &
689) geomechanics_regression_variable
690)
691) if (.not.associated(geomechanics_regression_variable)) return
692)
693) call GeomechanicsRegressionVariableDestroy( &
694) geomechanics_regression_variable%next)
695)
696) deallocate(geomechanics_regression_variable)
697) nullify(geomechanics_regression_variable)
698)
699) end subroutine GeomechanicsRegressionVariableDestroy
700)
701) ! ************************************************************************** !
702)
703) subroutine GeomechanicsRegressionDestroy(geomechanics_regression)
704) !
705) ! Destroys a geomechanics regression object
706) !
707) ! Author: Satish Karra
708) ! Date: 06/22/2016
709) !
710)
711) use Utility_module
712)
713) implicit none
714)
715) type(geomechanics_regression_type), pointer :: geomechanics_regression
716)
717) PetscErrorCode :: ierr
718)
719) if (.not.associated(geomechanics_regression)) return
720)
721) call GeomechanicsRegressionVariableDestroy( &
722) geomechanics_regression%variable_list)
723) call DeallocateArray(geomechanics_regression%natural_vertex_ids)
724) geomechanics_regression%num_vertices_per_process = 0
725) call DeallocateArray( &
726) geomechanics_regression%vertices_per_process_natural_ids)
727) if (geomechanics_regression%natural_vertex_id_vec /= 0) then
728) call VecDestroy(geomechanics_regression%natural_vertex_id_vec, &
729) ierr);CHKERRQ(ierr)
730) endif
731) if (geomechanics_regression%vertices_per_process_vec /= 0) then
732) call VecDestroy(geomechanics_regression%vertices_per_process_vec, &
733) ierr);CHKERRQ(ierr)
734) endif
735) if (geomechanics_regression%scatter_natural_vertex_id_gtos /= 0) then
736) call VecScatterDestroy( &
737) geomechanics_regression%scatter_natural_vertex_id_gtos, &
738) ierr);CHKERRQ(ierr)
739) endif
740) if (geomechanics_regression%scatter_vertices_per_process_gtos /= 0) then
741) call VecScatterDestroy( &
742) geomechanics_regression%scatter_vertices_per_process_gtos, &
743) ierr);CHKERRQ(ierr)
744) endif
745)
746) deallocate(geomechanics_regression)
747) nullify(geomechanics_regression)
748)
749) end subroutine GeomechanicsRegressionDestroy
750)
751) end module Geomechanics_Regression_module