/[MITgcm]/MITgcm_contrib/dgoldberg/streamice/streamice_cg_solve.F
ViewVC logotype

Diff of /MITgcm_contrib/dgoldberg/streamice/streamice_cg_solve.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.5 by dgoldberg, Thu Jul 26 16:13:18 2012 UTC revision 1.10 by dgoldberg, Wed Aug 27 19:29:13 2014 UTC
# Line 17  CBOP Line 17  CBOP
17       I                               A_vv,       ! section of matrix that multiplies v and projects on v       I                               A_vv,       ! section of matrix that multiplies v and projects on v
18       I                               tolerance,       I                               tolerance,
19       O                               iters,       O                               iters,
20         I                               maxIter,
21       I                               myThid )       I                               myThid )
22  C     /============================================================\  C     /============================================================\
23  C     | SUBROUTINE                                                 |    C     | SUBROUTINE                                                 |  
# Line 26  C     | Line 27  C     |
27  C     \============================================================/  C     \============================================================/
28        IMPLICIT NONE        IMPLICIT NONE
29    
 C     === Global variables ===  
30  #include "SIZE.h"  #include "SIZE.h"
31  #include "EEPARAMS.h"  #include "EEPARAMS.h"
32  #include "PARAMS.h"  #include "PARAMS.h"
33  #include "STREAMICE.h"  #include "STREAMICE.h"
34  #include "STREAMICE_CG.h"  #include "STREAMICE_CG.h"
35    
36    
37    
38    !#ifdef ALLOW_PETSC
39    !#include "finclude/petsc.h"
40    ! UNCOMMENT IF V3.0
41    !#include "finclude/petscvec.h"
42    !#include "finclude/petscmat.h"
43    !#include "finclude/petscksp.h"
44    !#include "finclude/petscpc.h"
45    !#endif
46    C     === Global variables ===
47    
48                
49  C     !INPUT/OUTPUT ARGUMENTS  C     !INPUT/OUTPUT ARGUMENTS
50  C     cg_Uin, cg_Vin - input and output velocities  C     cg_Uin, cg_Vin - input and output velocities
51  C     cg_Bu, cg_Bv - driving stress  C     cg_Bu, cg_Bv - driving stress
52        INTEGER myThid        INTEGER myThid
53        INTEGER iters        INTEGER iters
54          INTEGER maxIter
55        _RL tolerance        _RL tolerance
56        _RL cg_Uin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)        _RL cg_Uin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
57        _RL cg_Vin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)        _RL cg_Vin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
# Line 58  C     LOCAL VARIABLES Line 71  C     LOCAL VARIABLES
71        _RL dot_p2_tile (nSx,nSy)        _RL dot_p2_tile (nSx,nSy)
72        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
73    
74        iters = streamice_max_cg_iter  
75    !#ifdef ALLOW_PETSC
76    !      INTEGER indices(2*(snx*nsx*sny*nsy))
77    !      INTEGER n_dofs_cum_sum (0:nPx*nPy-1), idx(1)
78    !      _RL rhs_values(2*(snx*nsx*sny*nsy))
79    !      _RL solution_values(2*(snx*nsx*sny*nsy))
80    !      _RL mat_values (2*Nx*Ny,2*(snx*nsx*sny*nsy))
81    !      _RL mat_values (18,1), mat_val_return(1)
82    !      INTEGER indices_col(18)
83    !      INTEGER local_dofs, global_dofs, dof_index, dof_index_col
84    !      INTEGER local_offset
85    !      Mat matrix
86    !      KSP ksp
87    !      PC  pc
88    !      Vec rhs
89    !      Vec solution
90    !      PetscErrorCode ierr
91    !#ifdef ALLOW_USE_MPI
92    !      integer mpiRC, mpiMyWid
93    !#endif
94    !#endif
95    
96    
97  #ifdef ALLOW_STREAMICE  #ifdef ALLOW_STREAMICE
98    
99    
100    
101          CALL TIMER_START ('STREAMICE_CG_SOLVE',myThid)
102    #ifndef STREAMICE_SERIAL_TRISOLVE
103    
104    #ifdef ALLOW_PETSC
105    
106          CALL STREAMICE_CG_SOLVE_PETSC(
107         U         cg_Uin,     ! x-velocities
108         U         cg_Vin,     ! y-velocities
109         I         cg_Bu,      ! force in x dir
110         I         cg_Bv,      ! force in y dir
111         I         A_uu,       ! section of matrix that multiplies u and projects on u
112         I         A_uv,       ! section of matrix that multiplies v and projects on u
113         I         A_vu,       ! section of matrix that multiplies u and projects on v
114         I         A_vv,       ! section of matrix that multiplies v and projects on v
115         I         tolerance,
116         I         maxIter,
117         O         iters,
118         I         myThid )
119    
120    
121    #else  /* ALLOW_PETSC */
122    
123    
124          iters = maxIter
125        conv_flag = 0        conv_flag = 0
126    
127    
# Line 134  C     FIND INITIAL RESIDUAL, and initial Line 193  C     FIND INITIAL RESIDUAL, and initial
193         ENDDO         ENDDO
194        ENDDO        ENDDO
195    
         
   
196        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
197         DO bi = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
198          DO j=1,sNy          DO j=1,sNy
# Line 149  C     FIND INITIAL RESIDUAL, and initial Line 206  C     FIND INITIAL RESIDUAL, and initial
206         ENDDO         ENDDO
207        ENDDO        ENDDO
208    
209    
210        CALL GLOBAL_SUM_TILE_RL( dot_p1_tile, dot_p1, myThid )        CALL GLOBAL_SUM_TILE_RL( dot_p1_tile, dot_p1, myThid )
211        resid_0 = sqrt(dot_p1)        resid_0 = sqrt(dot_p1)
212    
213          DO bj = myByLo(myThid), myByHi(myThid)
214           DO bi = myBxLo(myThid), myBxHi(myThid)
215    
216           WRITE(msgBuf,'(A,I1,I1,E14.7)') 'CONJ GRAD INIT RESID LOCAL, ',
217         &         bi,bj, dot_p1_tile(bi,bj)
218           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
219         &                     SQUEEZE_RIGHT , 1)
220    
221           enddo
222          enddo
223    
224        WRITE(msgBuf,'(A,E14.7)') 'CONJ GRAD INIT RESID, ',        WRITE(msgBuf,'(A,E14.7)') 'CONJ GRAD INIT RESID, ',
225       &         resid_0       &         resid_0
226         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
# Line 200  c  !!!!!!!!!!!!!!!!!! Line 269  c  !!!!!!!!!!!!!!!!!!
269        
270  c  ! initially, b-grid data is valid up to 3 halo nodes out -- right? (check for MITgcm!!)  c  ! initially, b-grid data is valid up to 3 halo nodes out -- right? (check for MITgcm!!)
271    
272        print *, "BEGINNING MAIN CG LOOP"        WRITE(msgBuf,'(A)') 'BEGINNING MAIN CG LOOP'
273           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
274         &                     SQUEEZE_RIGHT , 1)
275    
276  !       IF(STREAMICE_construct_matrix) CALL STREAMICE_CG_MAKE_A(myThid)  !       IF(STREAMICE_construct_matrix) CALL STREAMICE_CG_MAKE_A(myThid)
277    
278    
279        do iter = 1, streamice_max_cg_iter        do iter = 1, maxIter
280         if (resid .gt. tolerance*resid_0) then         if (resid .gt. tolerance*resid_0) then
281    
282  c      to avoid using "exit"  c      to avoid using "exit"
# Line 429  c      to avoid using "exit" Line 500  c      to avoid using "exit"
500  c     to avoid using "exit"  c     to avoid using "exit"
501  c     if iters has reached max_iters there is no convergence  c     if iters has reached max_iters there is no convergence
502                
503        IF (iters .lt. streamice_max_cg_iter) THEN        IF (iters .lt. maxIter) THEN
504         conv_flag = 1         conv_flag = 1
505        ENDIF        ENDIF
506    
# Line 449  c     if iters has reached max_iters the Line 520  c     if iters has reached max_iters the
520  !       _EXCH_XY_RL( cg_Uin, myThid )  !       _EXCH_XY_RL( cg_Uin, myThid )
521  !       _EXCH_XY_RL( cg_Vin, myThid )          !       _EXCH_XY_RL( cg_Vin, myThid )        
522    
523  #endif  #endif /* ifndef ALLOW_PETSC */
       RETURN  
       END  
524    
525    #else /* STREAMICE_SERIAL_TRISOLVE */
526    
527                iters = 0
528    
529          CALL STREAMICE_TRIDIAG_SOLVE(
530         U                               cg_Uin,     ! x-velocities
531         U                               cg_Vin,
532         U                               cg_Bu,      ! force in x dir
533         I                               A_uu,       ! section of matrix that multiplies u and projects on u
534         I                               STREAMICE_umask,
535         I                               myThid )
536    
537    #endif
538    
539          CALL TIMER_STOP ('STREAMICE_CG_SOLVE',myThid)
540    
541    
542    #endif
543          RETURN
544          END
Legend:
Removed from v.1.5  
changed lines
  Added in v.1.10
  ViewVC Help
Powered by ViewVC 1.1.22