convergence.F90       coverage:  100.00 %func     21.38 %block


     1) module Convergence_module
     2) 
     3)   use Solver_module
     4)   use Option_module
     5)   use Grid_module
     6)   
     7)   use PFLOTRAN_Constants_module
     8) 
     9)   implicit none
    10) 
    11)   private
    12) 
    13) #include "petsc/finclude/petscsys.h"
    14) #include "petsc/finclude/petscvec.h"
    15) #include "petsc/finclude/petscvec.h90"
    16) #include "petsc/finclude/petscksp.h"
    17) #include "petsc/finclude/petscsnes.h"
    18) #include "petsc/finclude/petsclog.h"
    19) 
    20)   type, public :: convergence_context_type
    21)     type(solver_type), pointer :: solver
    22)     type(option_type), pointer :: option
    23)     type(grid_type), pointer :: grid
    24)   end type convergence_context_type
    25) 
    26) 
    27)   public :: ConvergenceContextCreate, ConvergenceTest, &
    28)             ConvergenceContextDestroy
    29)   
    30) contains
    31) 
    32) ! ************************************************************************** !
    33) 
    34) function ConvergenceContextCreate(solver,option,grid)
    35)   ! 
    36)   ! Creates a context containing pointer
    37)   ! for convergence subroutines
    38)   ! 
    39)   ! Author: Glenn Hammond
    40)   ! Date: 02/12/08
    41)   ! 
    42) 
    43)   implicit none
    44)   
    45)   type(convergence_context_type), pointer :: ConvergenceContextCreate
    46)   type(solver_type), pointer :: solver
    47)   type(option_type), pointer :: option
    48)   type(grid_type), pointer :: grid
    49)   
    50)   type(convergence_context_type), pointer :: context
    51)   
    52)   allocate(context)
    53)   context%solver => solver
    54)   context%option => option
    55)   context%grid => grid
    56) 
    57)   ConvergenceContextCreate => context
    58) 
    59) end function ConvergenceContextCreate
    60) 
    61) ! ************************************************************************** !
    62) 
    63) subroutine ConvergenceTest(snes_,i_iteration,xnorm,unorm,fnorm,reason,context, &
    64)                            ierr)
    65)   ! 
    66)   ! User defined convergence test
    67)   ! 
    68)   ! Author: Glenn Hammond
    69)   ! Date: 02/12/08
    70)   ! 
    71) 
    72)   implicit none
    73)   
    74)   SNES :: snes_
    75)   PetscInt :: i_iteration
    76)   PetscReal :: xnorm ! 2-norm of updated solution
    77)   PetscReal :: unorm ! 2-norm of update. PETSc refers to this as snorm
    78)   PetscReal :: fnorm ! 2-norm of updated residual
    79)   SNESConvergedReason :: reason
    80)   type(convergence_context_type) :: context
    81)   PetscErrorCode :: ierr
    82)   
    83)   type(solver_type), pointer :: solver
    84)   type(option_type), pointer :: option
    85)   type(grid_type), pointer :: grid
    86)   
    87)   Vec :: solution_vec
    88)   Vec :: update_vec
    89)   Vec :: residual_vec
    90)   PetscReal :: inorm_solution  !infinity norms
    91)   PetscReal :: inorm_update  
    92)   PetscReal :: inorm_residual  
    93)   
    94)   PetscInt :: i, ndof, max_index, min_index
    95)   PetscReal, allocatable :: fnorm_solution_stride(:)
    96)   PetscReal, allocatable :: fnorm_update_stride(:)
    97)   PetscReal, allocatable :: fnorm_residual_stride(:)
    98)   PetscReal, allocatable :: inorm_solution_stride(:)
    99)   PetscReal, allocatable :: inorm_update_stride(:)
   100)   PetscReal, allocatable :: inorm_residual_stride(:)
   101)   
   102)   PetscReal :: norm1_solution
   103)   PetscReal :: norm1_update
   104)   PetscReal :: norm1_residual
   105)   PetscReal, allocatable :: norm1_solution_stride(:)
   106)   PetscReal, allocatable :: norm1_update_stride(:)
   107)   PetscReal, allocatable :: norm1_residual_stride(:)
   108) 
   109)   KSP :: ksp
   110)   
   111)   PetscInt, allocatable :: imax_solution(:)
   112)   PetscInt, allocatable :: imax_update(:)
   113)   PetscInt, allocatable :: imax_residual(:)
   114)   PetscReal, allocatable :: max_solution_val(:)
   115)   PetscReal, allocatable :: max_update_val(:)
   116)   PetscReal, allocatable :: max_residual_val(:)
   117)   
   118)   PetscInt, allocatable :: imin_solution(:)
   119)   PetscInt, allocatable :: imin_update(:)
   120)   PetscInt, allocatable :: imin_residual(:)
   121)   PetscReal, allocatable :: min_solution_val(:)
   122)   PetscReal, allocatable :: min_update_val(:)
   123)   PetscReal, allocatable :: min_residual_val(:)
   124)   PetscReal, pointer :: vec_ptr(:)
   125)   
   126)   character(len=MAXSTRINGLENGTH) :: string, string2, string3, sec_string
   127)   PetscBool :: print_sol_norm_info = PETSC_FALSE
   128)   PetscBool :: print_upd_norm_info = PETSC_FALSE
   129)   PetscBool :: print_res_norm_info = PETSC_FALSE
   130)   PetscBool :: print_norm_by_dof_info = PETSC_FALSE
   131)   PetscBool :: print_max_val_and_loc_info = PETSC_FALSE
   132)   PetscBool :: print_1_norm_info = PETSC_FALSE
   133)   PetscBool :: print_2_norm_info = PETSC_FALSE
   134)   PetscBool :: print_inf_norm_info = PETSC_FALSE
   135) 
   136)   PetscInt :: sec_reason
   137) 
   138) !typedef enum {/* converged */
   139) !              SNES_CONVERGED_FNORM_ABS         =  2, /* F < F_minabs */
   140) !              SNES_CONVERGED_FNORM_RELATIVE    =  3, /* F < F_mintol*F_initial */
   141) !              SNES_CONVERGED_SNORM_RELATIVE    =  4, /* step size small */
   142) !              SNES_CONVERGED_ITS               =  5, /* maximum iterations reached */
   143) !              SNES_CONVERGED_TR_DELTA          =  7,
   144) !              /* diverged */
   145) !              SNES_DIVERGED_FUNCTION_DOMAIN    = -1,  
   146) !              SNES_DIVERGED_FUNCTION_COUNT     = -2,  
   147) !              SNES_DIVERGED_LINEAR_SOLVE       = -3, 
   148) !              SNES_DIVERGED_FNORM_NAN          = -4, 
   149) !              SNES_DIVERGED_MAX_IT             = -5,
   150) !              SNES_DIVERGED_LS_FAILURE         = -6,
   151) !              SNES_DIVERGED_LOCAL_MIN          = -8,  /* || J^T b || is small,
   152) !                                        implies converged to local minimum of F() */
   153) !              SNES_CONVERGED_ITERATING         =  0} SNESConvergedReason;
   154) 
   155)   solver => context%solver
   156)   option => context%option
   157)   grid => context%grid
   158) 
   159)   if (option%use_touch_options) then
   160)     string = 'detailed_convergence'
   161)     if (OptionCheckTouch(option,string)) then
   162)       if (solver%print_detailed_convergence) then
   163)         solver%print_detailed_convergence = PETSC_FALSE
   164)       else
   165)         solver%print_detailed_convergence = PETSC_TRUE
   166)       endif
   167)     endif
   168)   endif
   169) 
   170)   !geh: We must check the convergence here as i_iteration initializes
   171)   !     snes->ttol for subsequent iterations.
   172)   call SNESConvergedDefault(snes_,i_iteration,xnorm,unorm,fnorm,reason, &
   173)                             PETSC_NULL_OBJECT,ierr);CHKERRQ(ierr)
   174) #if 0
   175)   if (i_iteration == 0 .and. &
   176)       option%print_screen_flag .and. solver%print_convergence) then
   177)     write(*,'(i3," 2r:",es9.2)') i_iteration, fnorm
   178)   endif
   179) #endif
   180) 
   181)   ! for some reason (e.g. negative saturation/mole fraction in multiphase),
   182)   ! we are forcing extra newton iterations
   183)   if (option%force_newton_iteration) then
   184)     reason = 0
   185) !   reason = -1
   186)     return
   187)   endif
   188)   
   189) ! Checking if norm exceeds divergence tolerance
   190) !geh: inorm_residual is being used without being calculated.
   191) !      if (fnorm > solver%max_norm .or. unorm > solver%max_norm .or. &
   192) !        inorm_residual > solver%max_norm) then
   193)   
   194)   if (option%out_of_table) then
   195)     reason = -9
   196)   endif
   197)    
   198)   if (option%converged) then
   199)     reason = 12
   200)     ! set back to false
   201)     option%converged = PETSC_FALSE
   202)   endif
   203)     
   204) !  if (reason <= 0 .and. solver%check_infinity_norm) then
   205)   if (solver%check_infinity_norm) then
   206)   
   207)     call SNESGetFunction(snes_,residual_vec,PETSC_NULL_OBJECT, &
   208)                          PETSC_NULL_INTEGER,ierr);CHKERRQ(ierr)
   209) 
   210)     call VecNorm(residual_vec,NORM_INFINITY,inorm_residual,ierr);CHKERRQ(ierr)
   211) 
   212)     if (i_iteration > 0) then
   213)       call SNESGetSolutionUpdate(snes_,update_vec,ierr);CHKERRQ(ierr)
   214)       call VecNorm(update_vec,NORM_INFINITY,inorm_update,ierr);CHKERRQ(ierr)
   215)     else
   216)       inorm_update = 0.d0
   217)     endif
   218) 
   219)     if (inorm_residual < solver%newton_inf_res_tol) then
   220)       reason = 10
   221)     else
   222) !      if (reason > 0 .and. inorm_residual > 100.d0*solver%newton_inf_res_tol) &
   223) !        reason = 0
   224)     endif
   225) 
   226)     if (inorm_update < solver%newton_inf_upd_tol .and. i_iteration > 0) then
   227)       reason = 11
   228)     endif
   229)     
   230)     if (inorm_residual > solver%max_norm) then
   231)       reason = -10
   232)     endif
   233)  
   234)     ! This is to check if the secondary continuum residual convergences
   235)     ! for nonlinear problems specifically transport
   236)     if (solver%itype == TRANSPORT_CLASS .and. option%use_mc .and. &
   237)        reason > 0 .and. i_iteration > 0) then
   238)       if (option%infnorm_res_sec < solver%newton_inf_res_tol_sec) then
   239)         sec_reason = 1
   240)       else
   241)         reason = 0
   242)       endif
   243)     endif
   244)     
   245)     ! force the minimum number of iterations
   246)     if (i_iteration < solver%newton_min_iterations) then
   247)       reason = 0
   248)     endif
   249) 
   250)     if (option%print_screen_flag .and. solver%print_convergence) then
   251)       i = int(reason)
   252)       select case(i)
   253)         case(-10)
   254)           string = 'max_norm'
   255)         case(-9)
   256)           string = 'out_of_EOS_table'
   257)         case(2)
   258)           string = 'atol'
   259)         case(3)
   260)           string = 'rtol'
   261)         case(4)
   262)           string = 'stol'
   263)         case(10)
   264)           string = 'itol_res'
   265)         case(11)
   266)           string = 'itol_upd'
   267)         case(12)
   268)           string = 'itol_post_check'
   269)         case default
   270)           write(string,'(i3)') reason
   271)       end select
   272)       if (option%use_mc .and. option%ntrandof > 0 .and. solver%itype == &
   273)           TRANSPORT_CLASS) then
   274)         i = int(sec_reason) 
   275)         select case(i)
   276)           case(1)
   277)             sec_string = 'itol_res_sec'
   278)           case default
   279)             write(sec_string,'(i3)') sec_reason
   280)         end select
   281)         write(*,'(i3," 2r:",es9.2, &
   282)                 & " 2x:",es9.2, &
   283)                 & " 2u:",es9.2, &
   284)                 & " ir:",es9.2, &
   285)                 & " iu:",es9.2, &
   286)                 & " irsec:",es9.2, &
   287)                 & " rsn: ",a, ", ",a)') &
   288)                 i_iteration, fnorm, xnorm, unorm, inorm_residual, &
   289)                 inorm_update, option%infnorm_res_sec, &
   290)                 trim(string), trim(sec_string)
   291)       else
   292)         write(*,'(i3," 2r:",es9.2, &
   293)                 & " 2x:",es9.2, &
   294)                 & " 2u:",es9.2, &
   295)                 & " ir:",es9.2, &
   296)                 & " iu:",es9.2, &
   297)                 & " rsn: ",a)') &
   298)                 i_iteration, fnorm, xnorm, unorm, inorm_residual, &
   299)                 inorm_update, trim(string)        
   300)       endif
   301)     endif
   302)   else
   303)   
   304)     ! This is to check if the secondary continuum residual convergences
   305)     ! for nonlinear problems specifically transport
   306)     if (solver%itype == TRANSPORT_CLASS .and. option%use_mc .and. &
   307)        reason > 0 .and. i_iteration > 0) then
   308)       if (option%infnorm_res_sec < solver%newton_inf_res_tol_sec) then
   309)         reason = 13
   310)       else
   311)         reason = 0
   312)       endif
   313)     endif
   314)     
   315)     ! force the minimum number of iterations
   316)     if (i_iteration < solver%newton_min_iterations) then
   317)       reason = 0
   318)     endif
   319) 
   320)     if (option%print_screen_flag .and. solver%print_convergence) then
   321)       i = int(reason)
   322)       select case(i)
   323)         case(-9)
   324)           string = 'out_of_EOS_table'
   325)         case(2)
   326)           string = 'atol'
   327)         case(3)
   328)           string = 'rtol'
   329)         case(4)
   330)           string = 'stol'
   331)         case(10)
   332)           string = 'itol_res'
   333)         case(11)
   334)           string = 'itol_upd'
   335)         case(12)
   336)           string = 'itol_post_check'
   337)         case(13)
   338)           string = 'itol_res_sec'
   339)         case default
   340)           write(string,'(i3)') reason
   341)       end select
   342)       write(*,'(i3," 2r:",es10.2, &
   343)               & " 2u:",es10.2, &
   344)               & 32x, &
   345)               & " rsn: ",a)') i_iteration, fnorm, unorm, trim(string)
   346)       if (solver%print_linear_iterations) then
   347)         call KSPGetIterationNumber(solver%ksp,i,ierr);CHKERRQ(ierr)
   348)         write(option%io_buffer,'("   Linear Solver Iterations: ",i6)') i
   349)         call printMsg(option)
   350)       endif
   351)     endif
   352)   endif    
   353) 
   354)   if (solver%print_detailed_convergence) then
   355) 
   356)     call SNESGetSolution(snes_,solution_vec,ierr);CHKERRQ(ierr)
   357)     ! the ctx object should really be PETSC_NULL_OBJECT.  A bug in petsc
   358)     call SNESGetFunction(snes_,residual_vec,PETSC_NULL_OBJECT, &
   359)                          PETSC_NULL_INTEGER, &
   360)                          ierr);CHKERRQ(ierr)
   361)     call SNESGetSolutionUpdate(snes_,update_vec,ierr);CHKERRQ(ierr)
   362)     
   363)     ! infinity norms
   364)     call VecNorm(solution_vec,NORM_INFINITY,inorm_solution,ierr);CHKERRQ(ierr)
   365)     call VecNorm(update_vec,NORM_INFINITY,inorm_update,ierr);CHKERRQ(ierr)
   366)     call VecNorm(residual_vec,NORM_INFINITY,inorm_residual,ierr);CHKERRQ(ierr)
   367) 
   368)     call VecNorm(solution_vec,NORM_1,norm1_solution,ierr);CHKERRQ(ierr)
   369)     call VecNorm(update_vec,NORM_1,norm1_update,ierr);CHKERRQ(ierr)
   370)     call VecNorm(residual_vec,NORM_1,norm1_residual,ierr);CHKERRQ(ierr)
   371)     
   372)     call VecGetBlockSize(solution_vec,ndof,ierr);CHKERRQ(ierr)
   373)     
   374)     allocate(fnorm_solution_stride(ndof))
   375)     allocate(fnorm_update_stride(ndof))
   376)     allocate(fnorm_residual_stride(ndof))
   377)     allocate(inorm_solution_stride(ndof))
   378)     allocate(inorm_update_stride(ndof))
   379)     allocate(inorm_residual_stride(ndof))
   380)     allocate(norm1_solution_stride(ndof))
   381)     allocate(norm1_update_stride(ndof))
   382)     allocate(norm1_residual_stride(ndof))
   383)     
   384)     allocate(imax_solution(ndof))
   385)     allocate(imax_update(ndof))
   386)     allocate(imax_residual(ndof))
   387)     allocate(max_solution_val(ndof))
   388)     allocate(max_update_val(ndof))
   389)     allocate(max_residual_val(ndof))
   390) 
   391)     allocate(imin_solution(ndof))
   392)     allocate(imin_update(ndof))
   393)     allocate(imin_residual(ndof))
   394)     allocate(min_solution_val(ndof))
   395)     allocate(min_update_val(ndof))
   396)     allocate(min_residual_val(ndof))
   397) 
   398)     call VecStrideNormAll(solution_vec,NORM_1,norm1_solution_stride, &
   399)                           ierr);CHKERRQ(ierr)
   400)     call VecStrideNormAll(update_vec,NORM_1,norm1_update_stride, &
   401)                           ierr);CHKERRQ(ierr)
   402)     call VecStrideNormAll(residual_vec,NORM_1,norm1_residual_stride, &
   403)                           ierr);CHKERRQ(ierr)
   404)     call VecStrideNormAll(solution_vec,NORM_2,fnorm_solution_stride, &
   405)                           ierr);CHKERRQ(ierr)
   406)     call VecStrideNormAll(update_vec,NORM_2,fnorm_update_stride, &
   407)                           ierr);CHKERRQ(ierr)
   408)     call VecStrideNormAll(residual_vec,NORM_2,fnorm_residual_stride, &
   409)                           ierr);CHKERRQ(ierr)
   410)     call VecStrideNormAll(solution_vec,NORM_INFINITY,inorm_solution_stride, &
   411)                           ierr);CHKERRQ(ierr)
   412)     call VecStrideNormAll(update_vec,NORM_INFINITY,inorm_update_stride, &
   413)                           ierr);CHKERRQ(ierr)
   414)     call VecStrideNormAll(residual_vec,NORM_INFINITY,inorm_residual_stride, &
   415)                           ierr);CHKERRQ(ierr)
   416)     
   417)     ! can't use VecStrideMaxAll since the index location is not currently supported.
   418)     do i=1,ndof
   419)       call VecStrideMax(solution_vec,i-1,imax_solution(i),max_solution_val(i), &
   420)                         ierr);CHKERRQ(ierr)
   421)       call VecStrideMax(update_vec,i-1,imax_update(i),max_update_val(i), &
   422)                         ierr);CHKERRQ(ierr)
   423)       call VecStrideMax(residual_vec,i-1,imax_residual(i),max_residual_val(i), &
   424)                         ierr);CHKERRQ(ierr)
   425)       ! tweak the index to get the cell id from the mdof vector
   426)       imax_solution(i) = GridIndexToCellID(solution_vec,imax_solution(i),grid,GLOBAL)
   427)       imax_update(i) = GridIndexToCellID(update_vec,imax_update(i),grid,GLOBAL)
   428)       imax_residual(i) = GridIndexToCellID(residual_vec,imax_residual(i),grid,GLOBAL)
   429) !      imax_solution(i) = imax_solution(i)/ndof
   430) !      imax_update(i) = imax_update(i)/ndof
   431) !      imax_residual(i) = imax_residual(i)/ndof
   432)     enddo
   433) 
   434)     do i=1,ndof
   435)       call VecStrideMin(solution_vec,i-1,imin_solution(i),min_solution_val(i), &
   436)                         ierr);CHKERRQ(ierr)
   437)       call VecStrideMin(update_vec,i-1,imin_update(i),min_update_val(i), &
   438)                         ierr);CHKERRQ(ierr)
   439)       call VecStrideMin(residual_vec,i-1,imin_residual(i),min_residual_val(i), &
   440)                         ierr);CHKERRQ(ierr)
   441)       ! tweak the index to get the cell id from the mdof vector
   442)       imin_solution(i) = GridIndexToCellID(solution_vec,imin_solution(i),grid,GLOBAL)
   443)       imin_update(i) = GridIndexToCellID(update_vec,imax_update(i),grid,GLOBAL)
   444)       imin_residual(i) = GridIndexToCellID(residual_vec,imin_residual(i),grid,GLOBAL)
   445) !      imin_solution(i) = imin_solution(i)/ndof
   446) !      imin_update(i) = imin_update(i)/ndof
   447) !      imin_residual(i) = imin_residual(i)/ndof
   448)     enddo
   449) 
   450)     if (option%print_screen_flag) then
   451)       select case(reason)
   452)         case (10)
   453)           string = "CONVERGED_USER_NORM_INF_REL"
   454)         case (11)
   455)           string = "CONVERGED_USER_NORM_INF_UPD"
   456)         case(SNES_CONVERGED_FNORM_ABS)
   457)           string = "SNES_CONVERGED_FNORM_ABS"
   458)         case(SNES_CONVERGED_FNORM_RELATIVE)
   459)           string = "SNES_CONVERGED_FNORM_RELATIVE"
   460)         case(SNES_CONVERGED_SNORM_RELATIVE)
   461)           string = "SNES_CONVERGED_SNORM_RELATIVE"
   462)         case(SNES_CONVERGED_ITS)
   463)           string = "SNES_CONVERGED_ITS"
   464)         case(SNES_CONVERGED_TR_DELTA)
   465)           string = "SNES_CONVERGED_TR_DELTA"
   466)   !      case(SNES_DIVERGED_FUNCTION_DOMAIN)
   467)   !        string = "SNES_DIVERGED_FUNCTION_DOMAIN"
   468)         case(SNES_DIVERGED_FUNCTION_COUNT)
   469)           string = "SNES_DIVERGED_FUNCTION_COUNT"
   470)         case(SNES_DIVERGED_LINEAR_SOLVE)
   471)           string = "SNES_DIVERGED_LINEAR_SOLVE"
   472)         case(SNES_DIVERGED_FNORM_NAN)
   473)           string = "SNES_DIVERGED_FNORM_NAN"
   474)         case(SNES_DIVERGED_MAX_IT)
   475)           string = "SNES_DIVERGED_MAX_IT"
   476)         case(SNES_DIVERGED_LINE_SEARCH)
   477)           string = "SNES_DIVERGED_LINE_SEARCH"
   478)         case(SNES_DIVERGED_LOCAL_MIN)
   479)           string = "SNES_DIVERGED_LOCAL_MIN"
   480)         case(SNES_CONVERGED_ITERATING)
   481)           string = "SNES_CONVERGED_ITERATING"
   482)         case default
   483)           string = "UNKNOWN"
   484)       end select
   485) 
   486)       ! uncomment the lines below to determine data printed
   487)       
   488)       print_sol_norm_info = PETSC_TRUE  ! solution_vec norm information
   489)       print_upd_norm_info = PETSC_TRUE  ! update_vec norm information
   490)       print_res_norm_info = PETSC_TRUE  ! residual_vec norm information
   491)     
   492)       !print_norm_by_dof_info = PETSC_TRUE
   493)       print_max_val_and_loc_info = PETSC_TRUE
   494) 
   495)       !print_1_norm_info = PETSC_TRUE
   496)       print_2_norm_info = PETSC_TRUE
   497)       print_inf_norm_info = PETSC_TRUE
   498) 
   499)       print *
   500)       print *, 'reason: ', reason, ' - ', trim(string)
   501)       print *, 'SNES iteration :', i_iteration
   502)       call SNESGetKSP(snes_,ksp,ierr);CHKERRQ(ierr)
   503)       call KSPGetIterationNumber(ksp,i,ierr);CHKERRQ(ierr)
   504)       print *, 'KSP iterations :', i
   505)       if (print_1_norm_info) then
   506)         if (print_sol_norm_info) print *, 'norm_1_solution:   ', norm1_solution
   507)         if (print_upd_norm_info) print *, 'norm_1_update:     ', norm1_update
   508)         if (print_res_norm_info) print *, 'norm_1_residual:   ', norm1_residual
   509)       endif
   510)       if (print_2_norm_info) then
   511)         if (print_sol_norm_info) print *, 'norm_2_solution:   ', fnorm_solution_stride
   512)         if (print_upd_norm_info) print *, 'norm_2_update:     ', fnorm_update_stride
   513)         if (print_res_norm_info) print *, 'norm_2_residual:   ', fnorm_residual_stride
   514)       endif
   515)       if (print_inf_norm_info) then
   516)         if (print_sol_norm_info) print *, 'norm_inf_solution: ', inorm_solution
   517)         if (print_upd_norm_info) print *, 'norm_inf_update:   ', inorm_update
   518)         if (print_res_norm_info) print *, 'norm_inf_residual: ', inorm_residual
   519)       endif
   520)       if (print_max_val_and_loc_info) then
   521)         print *, 'max/min locations (zero-based index) by dof:'
   522)         do i=1,ndof
   523)           print *, '  dof: ', i
   524)           if (print_sol_norm_info) then
   525)             print *, '    solution_vec max: ', imax_solution(i), max_solution_val(i)
   526)             print *, '    solution_vec min: ', imin_solution(i), min_solution_val(i)
   527)           endif
   528)           if (print_upd_norm_info) then ! since update is -dx, need to invert
   529)             print *, '    update_vec max:   ', imin_update(i), -min_update_val(i)
   530)             print *, '    update_vec min:   ', imax_update(i), -max_update_val(i)
   531)           endif
   532)           if (print_res_norm_info) then
   533)             print *, '    residual_vec max: ', imax_residual(i), max_residual_val(i)
   534)             print *, '    residual_vec min: ', imin_residual(i), min_residual_val(i)
   535)           endif
   536)         enddo
   537)       endif
   538)       if (print_norm_by_dof_info) then
   539)         print *, 'norm by dof:'
   540)         do i=1,ndof
   541)           print *, '  dof: ', i
   542)           if (print_sol_norm_info) then
   543)             if (print_1_norm_info) &
   544)               print *, '    norm_1_solution:   ', norm1_solution_stride(i)
   545)             if (print_2_norm_info) &
   546)               print *, '    norm_2_solution:   ', fnorm_solution_stride(i)
   547)             if (print_inf_norm_info) &
   548)               print *, '    norm_inf_solution: ', inorm_solution_stride(i)
   549)             if (print_1_norm_info .or. print_2_norm_info .or. &
   550)                 print_inf_norm_info) print *, '    -'
   551)           endif
   552)           if (print_upd_norm_info) then
   553)             if (print_1_norm_info) &
   554)               print *, '    norm_1_update:   ', norm1_update_stride(i)
   555)             if (print_2_norm_info) &
   556)               print *, '    norm_2_update:   ', fnorm_update_stride(i)
   557)             if (print_inf_norm_info) &
   558)               print *, '    norm_inf_update: ', inorm_update_stride(i)
   559)             if (print_1_norm_info .or. print_2_norm_info .or. &
   560)                 print_inf_norm_info) print *, '    -'
   561)           endif
   562)           if (print_res_norm_info) then
   563)             if (print_1_norm_info) &
   564)               print *, '    norm_1_residual:   ', norm1_residual_stride(i)
   565)             if (print_2_norm_info) &
   566)               print *, '    norm_2_residual:   ', fnorm_residual_stride(i)
   567)             if (print_inf_norm_info) &
   568)               print *, '    norm_inf_residual: ', inorm_residual_stride(i)
   569)           endif
   570)         enddo
   571)       endif
   572)       print *
   573)     endif
   574)     
   575)     deallocate(fnorm_solution_stride)
   576)     deallocate(fnorm_update_stride)
   577)     deallocate(fnorm_residual_stride)
   578)     deallocate(inorm_solution_stride)
   579)     deallocate(inorm_update_stride)
   580)     deallocate(inorm_residual_stride)
   581)     deallocate(norm1_solution_stride)
   582)     deallocate(norm1_update_stride)
   583)     deallocate(norm1_residual_stride)
   584)     
   585)     deallocate(imax_solution)
   586)     deallocate(imax_update)
   587)     deallocate(imax_residual)
   588)     deallocate(max_solution_val)
   589)     deallocate(max_update_val)
   590)     deallocate(max_residual_val)
   591) 
   592)     deallocate(imin_solution)
   593)     deallocate(imin_update)
   594)     deallocate(imin_residual)
   595)     deallocate(min_solution_val)
   596)     deallocate(min_update_val)
   597)     deallocate(min_residual_val)
   598)     
   599)   endif
   600)   
   601) end subroutine ConvergenceTest
   602) 
   603) ! ************************************************************************** !
   604) 
   605) subroutine ConvergenceContextDestroy(context)
   606)   ! 
   607)   ! Destroy context
   608)   ! 
   609)   ! Author: Glenn Hammond
   610)   ! Date: 02/12/08
   611)   ! 
   612) 
   613)   implicit none
   614)   
   615)   type(convergence_context_type), pointer :: context
   616)   
   617)   if (.not.associated(context)) return
   618)   
   619)   nullify(context%solver)
   620)   nullify(context%option)
   621)   nullify(context%grid)
   622)   
   623)   deallocate(context)
   624)   nullify(context)
   625) 
   626) end subroutine ConvergenceContextDestroy
   627) 
   628) end module Convergence_module

generated by
Intel(R) C++/Fortran Compiler code-coverage tool
Web-Page Owner: Nobody