GlbThermCap.f90

Go to the documentation of this file.

!!****
!!
!!  NAME
!!     GlbThermCap
!!
!!  FUNCTION
!!     This subroutine constructs the global consistent capacitance matrix,
!!     then puts it into CRS format.  The matrix information is stored in
!!     global variables: rp_c, aval_c, cval_c, nnz_c, nrows_ktc, nstart_ktc
!!
!!  INPUTS
!!     none
!!
!!  OUTPUTS
!!     none
!!
!!  USES
!!     none
!!
!!****

SUBROUTINE glbthermcap(global)

  USE implicit_global
  USE comp_row_global
  USE rocstar_rocfraccomm

  IMPLICIT NONE

  include 'mpif.h'

  TYPE(rocfrac_global) :: global

  INTEGER :: i, j, m, p, inode, jnode, counter
  
  REAL(kind=wp) :: tempcval
  INTEGER, ALLOCATABLE, DIMENSION(:) :: amounttoreceive, amounttosend
  
  REAL(kind=wp), ALLOCATABLE, DIMENSION(:) :: bufsnd
  
  ! ... Gauss point coordinates
  REAL(kind=wp), DIMENSION(1:3,1:8) :: ri = reshape( &
       (/-0.577350269189626,-0.577350269189626,-0.577350269189626, &
          0.577350269189626,-0.577350269189626,-0.577350269189626, &
          0.577350269189626, 0.577350269189626,-0.577350269189626, &
         -0.577350269189626, 0.577350269189626,-0.577350269189626, &
         -0.577350269189626,-0.577350269189626, 0.577350269189626, &
          0.577350269189626,-0.577350269189626, 0.577350269189626, &
          0.577350269189626, 0.577350269189626, 0.577350269189626, &
         -0.577350269189626, 0.577350269189626, 0.577350269189626/),(/3,8/) )



  ! ... Construct the global capacitance matrix
  IF(myid == 0) print*,'CONSTRUCTING THE CAPACITANCE MATRIX'
  CALL mpi_barrier(rocstar_communicator,ierr)

  IF ( myid == 0) print*,'RocTherm, GlbThermCap :: before LocThermCap'

  ! ... Construct the capacitance matrix as a consistent matrix
  CALL locthermcap_v3d8(global%NumElVol,global%NumNp,global%NumMatVol,global%MeshCoor,global%ElConnVol,global%MatIdVol,ri,global%rho,global%Cp,global%ElConnVol)
  
  IF(myid == 0) print*,'FINISHED CONSTRUCTING THE CAPACITANCE MATRIX'
  CALL mpi_barrier(rocstar_communicator,ierr)
 

  ! ... number of non-zeros and CRS vectors ' '_temp computed by LocThermCap_v3d8
  nnz_c = nnz_temp
  ALLOCATE(cval_c(1:nnz_c))
  ALLOCATE(aval_c(1:nnz_c))
  ALLOCATE(rp_c(1:gnumnp+1))
  rp_c = rp_temp
  cval_c = cval_temp
  aval_c = aval_temp
  DEALLOCATE(rp_temp)
  DEALLOCATE(cval_temp)
  DEALLOCATE(aval_temp)


  IF(myid == 0) print*,'COMMUNICATING CAPACITANCE MATRIX TO ALL PROCESSES'

  ! ... Communicate mass to the other procs
  ! ... boundary nodes on this process are communicated to the processes
  ! ... to which they are assigned
  IF (nprocs > 1) THEN
     CALL mpi_barrier(rocstar_communicator,ierr)
     ! ... initialize request and status variables for MPI sends and receives
     ALLOCATE(req_rcv(1:nprocs))
     ALLOCATE(req_snd(1:nprocs))
     ALLOCATE(stat_snd(1:mpi_status_size,1:nprocs))
     ALLOCATE(stat_rcv(1:mpi_status_size,1:nprocs))
     DO i = 1, nprocs
        req_rcv(i) = 0
        req_snd(i) = 0
        DO j = 1, mpi_status_size
           stat_rcv(1,i) = 0
           stat_snd(j,i) = 0
        ENDDO
     ENDDO
     
     ! ... number of nodes each process will send to this process
     ! ... (nodes owned by this process but in the partition of other processes)
     ALLOCATE(amounttoreceive(1:numcommprocsfrom1))
     DO i = 1, numcommprocsfrom1
        CALL mpi_irecv(amounttoreceive(i),1, &
             mpi_integer,commprocsfrom1(i),10,rocstar_communicator, &
             req_rcv(i),ierr)
     ENDDO
print*,'RocTherm, GlbCap :: after MPI_IRECV',myid
     ! ... number of nodes each process will be sent from this process
     ! ... (nodes owned by other processes but in the partition of this one)
     ALLOCATE(amounttosend(1:numcommprocs1))
     DO i = 1, numcommprocs1
        counter = 0
        ! ... loop through the nodes going to be sent to process i
        DO j = 1, numcommnodes1(i)
           ! ... inode is the global node number corresponding to the jth node going to process i
           ! ... from this process
           inode = local2global(commnodes1(i,j))
           DO jnode = 1, gnumnp
              tempcval = 0.0
              ! ... return the value (tempCval) stored in CRS at explicit matrix
              ! ... location (inode, jnode)
              CALL comp_row_getval(gnumnp,gnumnp,nnz_c,1,rp_c,cval_c,aval_c, &
                   inode,jnode,tempcval)



              ! ... if there is a value at this node, then add 1 to the ammount to send
              ! ... to process i.
              IF (tempcval /= 0.0) THEN
                 counter = counter + 1
              ENDIF
           ENDDO
        ENDDO

        ! ... amount to send process i is the sum of all non-zeros in the loop above 
        amounttosend(i) = counter

        CALL mpi_send(amounttosend(i),1,mpi_integer, &
             commprocs1(i),10,rocstar_communicator,ierr)
     ENDDO

     CALL mpi_waitall(nprocs,req_rcv,stat_rcv,ierr)


     ! ... Communicate parts of the capcitance matrix to other procs

     ! ... receive information from other processes
     ! ... (nodes owned by this process but in the partition of other processes)
     CALL mpi_barrier(rocstar_communicator,ierr)
     ALLOCATE(frmproc(1:numcommprocsfrom1))
     DO i = 1, numcommprocsfrom1
        ALLOCATE(frmproc(i)%rcvbuf(1:3*amounttoreceive(i)))
        CALL mpi_irecv(frmproc(i)%rcvbuf(1),3*amounttoreceive(i), &
             mpi_double_precision,commprocsfrom1(i),10,rocstar_communicator, &
             req_rcv(i),ierr)
     ENDDO
     
     ! ... send information to other processes
     ! ... (nodes owned by other processes but in the partition of this one)    
     DO i = 1, numcommprocs1
        ALLOCATE(bufsnd(3*amounttosend(i)))
        bufsnd(:) = 0.0
        counter = 0
        ! ... loop through the nodes going to be sent to process i
        DO j = 1, numcommnodes1(i)
           ! ... inode is the global node number corresponding to the jth node going to process i
           ! ... from this process
           inode = local2global(commnodes1(i,j))
           DO jnode = 1, gnumnp
              tempcval = 0.0
              ! ... return the value (tempCval) stored in CRS at explicit matrix
              ! ... location (inode, jnode)
              CALL comp_row_getval(gnumnp,gnumnp,nnz_c,1,rp_c,cval_c,aval_c, &
                   inode,jnode,tempcval)

              ! ... if there is a value at this node, put it in the location '3*counter'
              ! ... (preceeded by the row and column locations (fortran)) in the buffer to
              ! ... be sent to process i.
              IF (tempcval /= 0.0) THEN
                 counter = counter + 1
                 bufsnd(3*counter-2) = inode
                 bufsnd(3*counter-1) = jnode
                 bufsnd(3*counter) = tempcval
              ENDIF
           ENDDO
        ENDDO

        CALL mpi_send(bufsnd,3*amounttosend(i),mpi_double_precision, &
             commprocs1(i),10,rocstar_communicator,ierr)

        DEALLOCATE(bufsnd)
     ENDDO

     CALL mpi_waitall(nprocs,req_rcv,stat_rcv,ierr)

     DEALLOCATE(req_rcv)
     DEALLOCATE(req_snd)
     DEALLOCATE(stat_snd)
     DEALLOCATE(stat_rcv)
  ENDIF

  IF(myid == 0) print*,'FINISHED COMMUNICATING CAPACITANCE MATRIX TO ALL PROCESSES'

  ! ... add values from other processes into this processes capacitance matrix
  IF (nprocs > 1) THEN
     DO i = 1, numcommprocsfrom1
        DO j = 1, amounttoreceive(i)
           CALL comp_row_addval(gnumnp,gnumnp,nnz_c,1,rp_c,cval_c,aval_c, &
                int(frmproc(i)%rcvbuf(3*j-2)),int(frmproc(i)%rcvbuf(3*j-1)),frmproc(i)%rcvbuf(3*j))

           ! ... if after the addition, the number of non-zeroes has changed
           ! ... reallocate the aval and cval vectors to accomodate
           IF (nnz_temp /= nnz_c) THEN
              DEALLOCATE(cval_c,aval_c)
              ALLOCATE(cval_c(nnz_temp),aval_c(nnz_temp))
              nnz_c = nnz_temp
              rp_c = rp_temp
              cval_c = cval_temp
           ENDIF
           aval_c = aval_temp
           DEALLOCATE(rp_temp,cval_temp,aval_temp)
        ENDDO
     ENDDO
     DEALLOCATE(frmproc)
     DEALLOCATE(amounttoreceive)
     DEALLOCATE(amounttosend)
  ENDIF

  ! ... Resize the matrix to the size needed on this proc
  CALL comp_row_resize(gnumnp,gnumnp,nnz_c,1,rp_c,cval_c,aval_c,nstart_ktc,nrows_ktc)
  DEALLOCATE(rp_c,cval_c,aval_c)
  ALLOCATE(rp_c(1:nrows_ktc+1))
  ALLOCATE(cval_c(1:nnz_temp))
  ALLOCATE(aval_c(1:nnz_temp))
  nnz_c = nnz_temp
  rp_c = rp_temp
  cval_c = cval_temp
  aval_c = aval_temp
  DEALLOCATE(rp_temp,cval_temp,aval_temp)
  

END SUBROUTINE glbthermcap