LCOV - code coverage report
Current view: top level - physics/clubb/src/CLUBB_core - matrix_solver_wrapper.F90 (source / functions) Hit Total Coverage
Test: coverage.info Lines: 40 95 42.1 %
Date: 2024-12-17 17:57:11 Functions: 4 5 80.0 %

          Line data    Source code
       1             : !-----------------------------------------------------------------------
       2             : ! $Id$
       3             : !===============================================================================
       4             : module matrix_solver_wrapper
       5             : 
       6             :   use clubb_precision, only: &
       7             :     core_rknd ! Variable(s)
       8             : 
       9             :   use error_code, only: &
      10             :     clubb_at_least_debug_level, & ! Procedure
      11             :     err_code,                   & ! Error indicator
      12             :     clubb_no_error,             & ! Constant
      13             :     clubb_fatal_error
      14             : 
      15             :   use constants_clubb, only: &
      16             :     fstderr     ! Constant(s)
      17             : 
      18             :   use model_flags, only: &
      19             :     lapack ,        & ! Variable(s)
      20             :     penta_lu ,      & 
      21             :     tridiag_lu,     &
      22             :     penta_bicgstab    
      23             : 
      24             :   implicit none
      25             :   
      26             :   public :: band_solve, &
      27             :             tridiag_solve
      28             : 
      29             :   interface band_solve
      30             :       module procedure band_solve_single_rhs_multiple_lhs
      31             :       module procedure band_solve_multiple_rhs_lhs
      32             :   end interface
      33             : 
      34             :   interface tridiag_solve
      35             :       module procedure tridiag_solve_single_rhs_lhs
      36             :       module procedure tridiag_solve_single_rhs_multiple_lhs
      37             :       module procedure tridiag_solve_multiple_rhs_lhs
      38             :   end interface
      39             :             
      40             :   contains
      41             : 
      42             : 
      43             :   !-------------------------------------------------------------------
      44             :   !                     Band Solvers Procedures
      45             :   !-------------------------------------------------------------------
      46             : 
      47     8935056 :   subroutine band_solve_single_rhs_multiple_lhs( &
      48             :                 solve_name, penta_solve_method, & ! Intent(in)
      49             :                 ngrdcol, nsup, nsub, ndim,      & ! Intent(in)
      50             :                 old_soln,                       & ! Intent(in)
      51     8935056 :                 lhs, rhs,                       & ! Intent(inout)
      52     8935056 :                 soln, rcond )                     ! Intent(out)
      53             : 
      54             :     use lapack_wrap, only:  & 
      55             :       lapack_band_solve,  & ! Procedure(s)
      56             :       lapack_band_solvex
      57             : 
      58             :     use penta_lu_solvers, only: &
      59             :       penta_lu_solve    ! Procedure(s)
      60             : 
      61             :     implicit none
      62             : 
      63             :     ! ----------------------- Input Variables -----------------------
      64             :     character(len=*), intent(in) :: &
      65             :       solve_name
      66             : 
      67             :     integer, intent(in) :: &
      68             :       penta_solve_method
      69             : 
      70             :     integer, intent(in) :: & 
      71             :       ngrdcol,  & ! Number of grid columns
      72             :       nsup,     & ! Number of superdiagonals
      73             :       nsub,     & ! Number of subdiagonals
      74             :       ndim        ! The order of the LHS Matrix, i.e. the # of linear equations
      75             : 
      76             :     real( kind = core_rknd ), dimension(ngrdcol,ndim), intent(in) ::  & 
      77             :       old_soln ! Old soln, used as an initial guess in the bicgstab method
      78             : 
      79             :     real( kind = core_rknd ), dimension(nsup+nsub+1,ngrdcol,ndim), intent(inout) ::  & 
      80             :       lhs ! Left hand side
      81             : 
      82             :     real( kind = core_rknd ), dimension(ngrdcol,ndim), intent(inout) ::  & 
      83             :       rhs ! Right hand side(s)
      84             : 
      85             :     ! ----------------------- Output Variables -----------------------
      86             :     real( kind = core_rknd ), dimension(ngrdcol,ndim), intent(out) :: &
      87             :       soln
      88             : 
      89             :     ! ----------------------- Optional Out -----------------------
      90             : 
      91             :     ! The estimate of the reciprocal condition number of matrix
      92             :     ! after equilibration (if done).
      93             :     real( kind = core_rknd ), optional, dimension(ngrdcol), intent(out) ::  & 
      94             :       rcond
      95             : 
      96             :     ! ----------------------- Local Variables -----------------------
      97             :     real( kind = core_rknd ), dimension(nsup+nsub+1,ngrdcol,ndim) ::  & 
      98    17870112 :       lhs_copy ! Copy of left hand side
      99             : 
     100             :     real( kind = core_rknd ), dimension(ngrdcol,ndim) ::  & 
     101    17870112 :       rhs_copy ! Copy of right hand side
     102             : 
     103             :     real( kind = core_rknd ), dimension(ngrdcol,ndim) :: &
     104    17870112 :       dummy_soln
     105             : 
     106             :     integer :: i
     107             : 
     108             :     ! ----------------------- Begin Code -----------------------
     109             : 
     110     8935056 :     if ( present(rcond) ) then
     111             : 
     112             :       !$acc update host( lhs, rhs )
     113             : 
     114             :       ! Lapack overwrites lhs and rhs, so we'll give it copies of them.
     115           0 :       lhs_copy = lhs
     116           0 :       rhs_copy = rhs
     117             : 
     118             :       ! Perform LU decomp and solve system (LAPACK with diagnostics)
     119             :       ! Using dummy_soln, since we only want this routine for diagnostics
     120             :       call lapack_band_solvex( "xm_wpxp", nsup, nsub, & ! Intent(in) 
     121             :                                ndim, 1, ngrdcol,      & ! Intent(in) 
     122             :                                lhs_copy, rhs_copy,    & ! Intent(inout)
     123           0 :                                dummy_soln, rcond )      ! Intent(out)
     124             : 
     125             :       !$acc update device( rcond )
     126             : 
     127             :     end if
     128             : 
     129             : 
     130     8935056 :     if ( penta_solve_method == lapack ) then
     131             : 
     132             :       !$acc update host( lhs, rhs )
     133             : 
     134             :       ! Perform LU decomp and solve system (LAPACK)
     135             :       call lapack_band_solve( "xm_wpxp", nsup, nsub,  & ! Intent(in) 
     136             :                               ndim, 1, ngrdcol,       & ! Intent(in) 
     137             :                               lhs, rhs,               & ! Intent(inout)
     138     8935056 :                               soln )                    ! Intent(out)
     139             : 
     140             :       !$acc update device( soln )
     141             : 
     142           0 :     else if ( penta_solve_method == penta_lu ) then 
     143             : 
     144             :       ! Solve the system with a penta-diagonal specific LU decomp
     145             :       call penta_lu_solve( ndim, ngrdcol,         & ! Intent(in)
     146             :                            lhs(:,:,:), rhs(:,:),  & ! Intent(in)
     147           0 :                            soln(:,:) )              ! Intent(out)
     148             : 
     149             :     else
     150             : 
     151             :       ! The solve method should match one of the above
     152           0 :       if ( clubb_at_least_debug_level( 0 ) ) then
     153           0 :         write(fstderr,*) "Error in band_solve_single_rhs_multiple_lhs: "
     154           0 :         write(fstderr,*) "  no case for penta_solve_method = ", penta_solve_method
     155           0 :         err_code = clubb_fatal_error
     156             :       end if
     157             : 
     158             :     end if
     159             : 
     160     8935056 :     return
     161             : 
     162             :   end subroutine band_solve_single_rhs_multiple_lhs
     163             : 
     164     8935056 :   subroutine band_solve_multiple_rhs_lhs( &
     165             :                 solve_name, penta_solve_method,   & ! Intent(in)
     166             :                 ngrdcol, nsup, nsub, ndim, nrhs,  & ! Intent(in)
     167             :                 old_soln,                         & ! Intent(in)
     168     8935056 :                 lhs, rhs,                         & ! Intent(inout)
     169     8935056 :                 soln, rcond )                       ! Intent(out)
     170             : 
     171             :     use lapack_wrap, only:  & 
     172             :       lapack_band_solve,  & ! Procedure(s)
     173             :       lapack_band_solvex
     174             : 
     175             :     use penta_lu_solvers, only: &
     176             :       penta_lu_solve    ! Procedure(s)
     177             : 
     178             :     implicit none
     179             : 
     180             :     ! ----------------------- Input Variables -----------------------
     181             :     character(len=*), intent(in) :: &
     182             :       solve_name
     183             : 
     184             :     integer, intent(in) :: &
     185             :       penta_solve_method
     186             : 
     187             :     integer, intent(in) :: & 
     188             :       ngrdcol,  & ! Number of grid columns
     189             :       nsup,     & ! Number of superdiagonals
     190             :       nsub,     & ! Number of subdiagonals
     191             :       ndim,     & ! The order of the LHS Matrix, i.e. the # of linear equations
     192             :       nrhs        ! Number of RHS's to back substitute for
     193             : 
     194             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs), intent(in) ::  & 
     195             :       old_soln ! Old soln, used as an initial guess in the bicgstab method
     196             : 
     197             :     real( kind = core_rknd ), dimension(nsup+nsub+1,ngrdcol,ndim), intent(inout) ::  & 
     198             :       lhs ! Left hand side
     199             : 
     200             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs), intent(inout) ::  & 
     201             :       rhs ! Right hand side(s)
     202             : 
     203             :     ! ----------------------- Output Variables -----------------------
     204             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs), intent(out) :: &
     205             :       soln
     206             : 
     207             :     ! ----------------------- Optional Out -----------------------
     208             : 
     209             :     ! The estimate of the reciprocal condition number of matrix
     210             :     ! after equilibration (if done).
     211             :     real( kind = core_rknd ), optional, dimension(ngrdcol), intent(out) ::  & 
     212             :       rcond
     213             : 
     214             :     ! ----------------------- Local Variables -----------------------
     215             :     real( kind = core_rknd ), dimension(nsup+nsub+1,ngrdcol,ndim) ::  & 
     216    17870112 :       lhs_copy ! Copy of left hand side
     217             : 
     218             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs) ::  & 
     219    17870112 :       rhs_copy ! Copy of right hand side
     220             : 
     221             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs) :: &
     222    17870112 :       dummy_soln
     223             : 
     224             :     integer :: i
     225             : 
     226             :     ! ----------------------- Begin Code -----------------------
     227             : 
     228     8935056 :     if ( present(rcond) ) then
     229             : 
     230             :       !$acc update host( lhs, rhs )
     231             : 
     232             :       ! Lapack overwrites lhs and rhs, so we'll give it copies of them.
     233           0 :       lhs_copy = lhs
     234           0 :       rhs_copy = rhs
     235             : 
     236             :       ! Perform LU decomp and solve system (LAPACK with diagnostics)
     237             :       ! Using dummy_soln, since we only want this routine for diagnostics
     238             :       call lapack_band_solvex( "xm_wpxp", nsup, nsub, & ! Intent(in) 
     239             :                                ndim, nrhs, ngrdcol,   & ! Intent(in) 
     240             :                                lhs_copy, rhs_copy,    & ! Intent(inout)
     241           0 :                                dummy_soln, rcond )      ! Intent(out)
     242             : 
     243             :       !$acc update device( rcond )
     244             : 
     245             :     end if
     246             : 
     247             : 
     248     8935056 :     if ( penta_solve_method == lapack ) then
     249             : 
     250             :       !$acc update host( lhs, rhs )
     251             : 
     252             :       ! Perform LU decomp and solve system (LAPACK)
     253             :       call lapack_band_solve( "xm_wpxp", nsup, nsub,  & ! Intent(in) 
     254             :                               ndim, nrhs, ngrdcol,    & ! Intent(in) 
     255             :                               lhs, rhs,               & ! Intent(inout)
     256     8935056 :                               soln )                    ! Intent(out)
     257             : 
     258             :       !$acc update device( soln )
     259             : 
     260           0 :     else if ( penta_solve_method == penta_lu ) then 
     261             : 
     262             :       ! Solve the system with a penta-diagonal specific LU decomp
     263             :       call penta_lu_solve( ndim, nrhs, ngrdcol,     & ! Intent(in)
     264             :                            lhs(:,:,:), rhs(:,:,:),  & ! Intent(in)
     265           0 :                            soln(:,:,:) )              ! Intent(out)
     266             : 
     267             :     else
     268             : 
     269             :       ! The solve method should match one of the above
     270           0 :       if ( clubb_at_least_debug_level( 0 ) ) then
     271           0 :         write(fstderr,*) "Error in band_solve_multiple_rhs_lhs: "
     272           0 :         write(fstderr,*) "  no case for penta_solve_method = ", penta_solve_method
     273           0 :         err_code = clubb_fatal_error
     274             :       end if
     275             : 
     276             :     end if
     277             : 
     278     8935056 :     return
     279             : 
     280             :   end subroutine band_solve_multiple_rhs_lhs
     281             : 
     282             : 
     283             : 
     284             :   !-------------------------------------------------------------------
     285             :   !                    Tridiag Solver Procedures
     286             :   !-------------------------------------------------------------------
     287             : 
     288           0 :   subroutine tridiag_solve_single_rhs_lhs( &
     289             :                 solve_name, tridiag_solve_method, & ! Intent(in)
     290             :                 ndim,                             & ! Intent(in)
     291           0 :                 lhs, rhs,                         & ! Intent(inout)
     292           0 :                 soln, rcond )                       ! Intent(out)
     293             : 
     294             :     use lapack_wrap, only:  & 
     295             :       lapack_tridiag_solve,  & ! Procedure(s)
     296             :       lapack_tridiag_solvex
     297             : 
     298             :     use tridiag_lu_solvers, only: &
     299             :       tridiag_lu_solve    ! Procedure(s)
     300             : 
     301             :     implicit none
     302             : 
     303             :     ! ----------------------- Input Variables -----------------------
     304             :     character(len=*), intent(in) :: &
     305             :       solve_name
     306             : 
     307             :     integer, intent(in) :: &
     308             :       tridiag_solve_method
     309             : 
     310             :     integer, intent(in) :: & 
     311             :       ndim        ! The order of the LHS Matrix, i.e. the # of linear equations
     312             : 
     313             :     real( kind = core_rknd ), dimension(3,ndim), intent(inout) ::  & 
     314             :       lhs ! Left hand side
     315             : 
     316             :     real( kind = core_rknd ), dimension(ndim), intent(inout) ::  & 
     317             :       rhs ! Right hand side(s)
     318             : 
     319             :     ! ----------------------- Output Variables -----------------------
     320             :     real( kind = core_rknd ), dimension(ndim), intent(out) :: &
     321             :       soln
     322             : 
     323             :     ! ----------------------- Optional Out -----------------------
     324             : 
     325             :     ! The estimate of the reciprocal condition number of matrix
     326             :     ! after equilibration (if done).
     327             :     real( kind = core_rknd ), dimension(1), optional, intent(out) ::  & 
     328             :       rcond
     329             : 
     330             :     ! ----------------------- Local Variables -----------------------
     331             :     real( kind = core_rknd ), dimension(3,ndim) ::  & 
     332           0 :       lhs_copy ! Copy of left hand side
     333             : 
     334             :     real( kind = core_rknd ), dimension(ndim) ::  & 
     335           0 :       rhs_copy ! Copy of right hand side
     336             : 
     337             :     real( kind = core_rknd ), dimension(ndim) :: &
     338           0 :       dummy_soln
     339             : 
     340             :     ! ----------------------- Begin Code -----------------------
     341             : 
     342           0 :     if ( present(rcond) ) then
     343             : 
     344             :       ! Lapack overwrites lhs and rhs, so we'll give it copies of them.
     345           0 :       lhs_copy = lhs
     346           0 :       rhs_copy = rhs
     347             : 
     348             :       ! Perform LU decomp and solve system (LAPACK with diagnostics)
     349             :       call lapack_tridiag_solvex( solve_name, ndim, 1, 1, & ! Intent(in) 
     350             :                                   lhs_copy, rhs_copy,     & ! Intent(inout)
     351           0 :                                   dummy_soln, rcond )       ! Intent(out)
     352             :     end if
     353             : 
     354             : 
     355           0 :     if ( tridiag_solve_method == lapack ) then
     356             : 
     357             :       ! Perform LU decomp and solve system (LAPACK)
     358             :       call lapack_tridiag_solve( solve_name, ndim, 1, 1, & ! Intent(in) 
     359             :                                  lhs, rhs,               & ! Intent(inout)
     360           0 :                                  soln )                   ! Intent(out)
     361             : 
     362           0 :     else if ( tridiag_solve_method == tridiag_lu ) then
     363             : 
     364             :       ! Solve the system with a tri-diagonal specific LU decomp
     365             :       call tridiag_lu_solve( ndim,      & ! Intent(in)
     366             :                              lhs, rhs,  & ! Intent(in)
     367           0 :                              soln )       ! Intent(out)
     368             : 
     369             :     else
     370             : 
     371             :       ! The solve method should match one of the above
     372           0 :       if ( clubb_at_least_debug_level( 0 ) ) then
     373           0 :         write(fstderr,*) "Error in tridiag_solve_single_rhs_lhs: "
     374           0 :         write(fstderr,*) "  no case for tridiag_solve_method = ", tridiag_solve_method
     375           0 :         err_code = clubb_fatal_error
     376             :       end if
     377             : 
     378             :     end if
     379             : 
     380           0 :     return
     381             : 
     382             :   end subroutine tridiag_solve_single_rhs_lhs
     383             : 
     384             : 
     385      405213 :   subroutine tridiag_solve_single_rhs_multiple_lhs( &
     386             :                 solve_name, tridiag_solve_method, & ! Intent(in)
     387             :                 ngrdcol, ndim,                    & ! Intent(in)
     388      405213 :                 lhs, rhs,                         & ! Intent(inout)
     389      405213 :                 soln, rcond )                       ! Intent(out)
     390             : 
     391             :     use lapack_wrap, only:  & 
     392             :       lapack_tridiag_solve,  & ! Procedure(s)
     393             :       lapack_tridiag_solvex
     394             : 
     395             :     use tridiag_lu_solvers, only: &
     396             :       tridiag_lu_solve    ! Procedure(s)
     397             : 
     398             :     implicit none
     399             : 
     400             :     ! ----------------------- Input Variables -----------------------
     401             :     character(len=*), intent(in) :: &
     402             :       solve_name
     403             : 
     404             :     integer, intent(in) :: &
     405             :       tridiag_solve_method
     406             : 
     407             :     integer, intent(in) :: & 
     408             :       ngrdcol,  & ! Number of grid columns
     409             :       ndim        ! The order of the LHS Matrix, i.e. the # of linear equations
     410             : 
     411             :     real( kind = core_rknd ), dimension(3,ngrdcol,ndim), intent(inout) ::  & 
     412             :       lhs ! Left hand side
     413             : 
     414             :     real( kind = core_rknd ), dimension(ngrdcol,ndim), intent(inout) ::  & 
     415             :       rhs ! Right hand side(s)
     416             : 
     417             :     ! ----------------------- Output Variables -----------------------
     418             :     real( kind = core_rknd ), dimension(ngrdcol,ndim), intent(out) :: &
     419             :       soln
     420             : 
     421             :     ! ----------------------- Optional Out -----------------------
     422             : 
     423             :     ! The estimate of the reciprocal condition number of matrix
     424             :     ! after equilibration (if done).
     425             :     real( kind = core_rknd ), optional, dimension(ngrdcol), intent(out) ::  & 
     426             :       rcond
     427             : 
     428             :     ! ----------------------- Local Variables -----------------------
     429             :     real( kind = core_rknd ), dimension(3,ngrdcol,ndim) ::  & 
     430      810426 :       lhs_copy ! Copy of left hand side
     431             : 
     432             :     real( kind = core_rknd ), dimension(ngrdcol,ndim) ::  & 
     433      810426 :       rhs_copy ! Copy of right hand side
     434             : 
     435             :     real( kind = core_rknd ), dimension(ngrdcol,ndim) :: &
     436      810426 :       dummy_soln
     437             : 
     438             :     integer :: i
     439             : 
     440             :     ! ----------------------- Begin Code -----------------------
     441             : 
     442      405213 :     if ( present(rcond) ) then
     443             : 
     444             :       !$acc update host( lhs, rhs )
     445             : 
     446             :       ! Lapack overwrites lhs and rhs, so we'll give it copies of them.
     447           0 :       lhs_copy = lhs
     448           0 :       rhs_copy = rhs
     449             : 
     450             :       ! Perform LU decomp and solve system (LAPACK with diagnostics)
     451             :       call lapack_tridiag_solvex( solve_name, ndim, 1, ngrdcol, & ! Intent(in) 
     452             :                                   lhs_copy, rhs_copy,           & ! Intent(inout)
     453           0 :                                   dummy_soln, rcond )             ! Intent(out)
     454             :       
     455             :       !$acc update device( rcond )
     456             : 
     457             :     end if
     458             : 
     459             : 
     460      405213 :     if ( tridiag_solve_method == lapack ) then
     461             : 
     462             :       !$acc update host( lhs, rhs )
     463             : 
     464             :       ! Perform LU decomp and solve system (LAPACK)
     465             :       call lapack_tridiag_solve( solve_name, ndim, 1, ngrdcol, & ! Intent(in) 
     466             :                                  lhs, rhs,                     & ! Intent(inout)
     467      405213 :                                  soln )                          ! Intent(out)
     468             : 
     469             :       !$acc update device( soln )
     470             : 
     471           0 :     else if ( tridiag_solve_method == tridiag_lu ) then
     472             : 
     473             :       ! Solve the system with a tri-diagonal specific LU decomp
     474             :       call tridiag_lu_solve( ndim, ngrdcol, & ! Intent(in)
     475             :                              lhs, rhs,      & ! Intent(in)
     476           0 :                              soln )           ! Intent(out)
     477             :     else
     478             : 
     479             :       ! The solve method should match one of the above
     480           0 :       if ( clubb_at_least_debug_level( 0 ) ) then
     481           0 :         write(fstderr,*) "Error in tridiag_solve_single_rhs_multiple_lhs: "
     482           0 :         write(fstderr,*) "  no case for tridiag_solve_method = ", tridiag_solve_method
     483           0 :         err_code = clubb_fatal_error
     484             :       end if
     485             : 
     486             :     end if
     487             : 
     488      405213 :     return
     489             : 
     490             :   end subroutine tridiag_solve_single_rhs_multiple_lhs
     491             : 
     492    44675280 :   subroutine tridiag_solve_multiple_rhs_lhs( &
     493             :                 solve_name, tridiag_solve_method, & ! Intent(in)
     494             :                 ngrdcol, ndim, nrhs,              & ! Intent(in)
     495    44675280 :                 lhs, rhs,                         & ! Intent(inout)
     496    44675280 :                 soln, rcond )                       ! Intent(out)
     497             : 
     498             :     use lapack_wrap, only:  & 
     499             :       lapack_tridiag_solve,  & ! Procedure(s)
     500             :       lapack_tridiag_solvex
     501             : 
     502             :     use tridiag_lu_solvers, only: &
     503             :       tridiag_lu_solve    ! Procedure(s)
     504             : 
     505             :     implicit none
     506             : 
     507             :     ! ----------------------- Input Variables -----------------------
     508             :     character(len=*), intent(in) :: &
     509             :       solve_name
     510             : 
     511             :     integer, intent(in) :: &
     512             :       tridiag_solve_method
     513             : 
     514             :     integer, intent(in) :: & 
     515             :       ngrdcol,  & ! Number of grid columns
     516             :       ndim,     & ! The order of the LHS Matrix, i.e. the # of linear equations
     517             :       nrhs        ! Number of RHS's to back substitute for
     518             : 
     519             :     real( kind = core_rknd ), dimension(3,ngrdcol,ndim), intent(inout) ::  & 
     520             :       lhs ! Left hand side
     521             : 
     522             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs), intent(inout) ::  & 
     523             :       rhs ! Right hand side(s)
     524             : 
     525             :     ! ----------------------- Output Variables -----------------------
     526             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs), intent(out) :: &
     527             :       soln
     528             : 
     529             :     ! ----------------------- Optional Out -----------------------
     530             : 
     531             :     ! The estimate of the reciprocal condition number of matrix
     532             :     ! after equilibration (if done).
     533             :     real( kind = core_rknd ), optional, dimension(ngrdcol), intent(out) ::  & 
     534             :       rcond
     535             : 
     536             :     ! ----------------------- Local Variables -----------------------
     537             :     real( kind = core_rknd ), dimension(3,ngrdcol,ndim) ::  & 
     538    89350560 :       lhs_copy ! Copy of left hand side
     539             : 
     540             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs) ::  & 
     541    89350560 :       rhs_copy ! Copy of right hand side
     542             : 
     543             :     real( kind = core_rknd ), dimension(ngrdcol,ndim,nrhs) :: &
     544    89350560 :       dummy_soln
     545             : 
     546             :     integer :: i
     547             : 
     548             :     ! ----------------------- Begin Code -----------------------
     549             : 
     550    44675280 :     if ( present(rcond) ) then
     551             : 
     552             :       !$acc update host( lhs, rhs )
     553             : 
     554             :       ! Lapack overwrites lhs and rhs, so we'll give it copies of them.
     555           0 :       lhs_copy = lhs
     556           0 :       rhs_copy = rhs
     557             : 
     558             :       ! Perform LU decomp and solve system (LAPACK with diagnostics)
     559             :       call lapack_tridiag_solvex( solve_name, ndim, nrhs, ngrdcol,  & ! Intent(in) 
     560             :                                   lhs_copy, rhs_copy,               & ! Intent(inout)
     561           0 :                                   dummy_soln, rcond )                 ! Intent(out)
     562             :       
     563             :       !$acc update device( rcond )
     564             : 
     565             :     end if
     566             : 
     567             : 
     568    44675280 :     if ( tridiag_solve_method == lapack ) then
     569             : 
     570             :       !$acc update host( lhs, rhs )
     571             : 
     572             :       ! Perform LU decomp and solve system (LAPACK)
     573             :       call lapack_tridiag_solve( solve_name, ndim, nrhs, ngrdcol, & ! Intent(in) 
     574             :                                  lhs, rhs,                        & ! Intent(inout)
     575    44675280 :                                  soln )                             ! Intent(out)
     576             : 
     577             :       !$acc update device( soln )
     578             : 
     579           0 :     else if ( tridiag_solve_method == tridiag_lu ) then
     580             : 
     581             :       ! Solve the system with a tri-diagonal specific LU decomp
     582             :       call tridiag_lu_solve( ndim, nrhs, ngrdcol, & ! Intent(in)
     583             :                              lhs, rhs,            & ! Intent(in)
     584           0 :                              soln )                 ! Intent(out)
     585             :     else
     586             : 
     587             :       ! The solve method should match one of the above
     588           0 :       if ( clubb_at_least_debug_level( 0 ) ) then
     589           0 :         write(fstderr,*) "Error in tridiag_solve_multiple_rhs_lhs: "
     590           0 :         write(fstderr,*) "  no case for tridiag_solve_method = ", tridiag_solve_method
     591           0 :         err_code = clubb_fatal_error
     592             :       end if
     593             : 
     594             :     end if
     595             : 
     596    44675280 :     return
     597             : 
     598             :   end subroutine tridiag_solve_multiple_rhs_lhs
     599             : 
     600             : 
     601             : end module matrix_solver_wrapper

Generated by: LCOV version 1.14