Commit e8b39720 authored by Serhiy Mochalskyy's avatar Serhiy Mochalskyy
Browse files

corrected output subroutine

parent 63e29e00
Pipeline #4484 skipped
......@@ -38,7 +38,7 @@ if(INFO_B .NE. 0) then
endif
CALL DESCINIT(DESCC, n_w,nd_bez, NB, NB, 0, 0, CONTEXT, LDA_ye, INFO_C )
CALL DESCINIT(DESC_ye, n_w,nd_bez, NB, NB, 0, 0, CONTEXT, LDA_ye, INFO_C )
if(INFO_C .NE. 0) then
write(6,*) "Something is wrong in a_ye_computing CALL DESCINIT DESCC, INFO_C=",INFO_C
stop
......@@ -46,7 +46,7 @@ endif
CALL PDGEMM('T','N', n_w ,nd_bez, n_w , 1.0D0 , s_ww_loc ,1,1,DESCA, a_we_loc,1,1,DESCB,&
1.0D0, a_ye_loc,1 , 1 , DESCC)
1.0D0, a_ye_loc,1 , 1 , DESC_ye)
do i=1,n_w
......
......@@ -80,7 +80,7 @@ subroutine input
call MPI_BCAST(NB, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ier)
call MPI_BCAST(ORFAC, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
call MPI_BCAST(lwork_cooficient, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ier)
!===================================================================================
!==================================================================================
if(rank==0) then
write(outp,*) ' ScaLAPACK Input parameters:'
write(outp,'(A,I4)') ' NB =', NB
......@@ -88,14 +88,16 @@ subroutine input
write(outp,'(A,I4)')
endif
!Only one task read file type because only one tast will write afterwards data to the file.
if(rank==0) read(inp, params_output) ! Read namelist from STDIN
!==================================================================================
call MPI_BCAST(format_type, len(format_type), MPI_CHARACTER,0, MPI_COMM_WORLD, ier)
!==================================================================================
if(rank==0) then
write(outp,'(A,I4)')
write(outp,*) 'Output starwall-response.dat:'
write(outp,'(A,A)') ' file format = ',format_type
endif
! --- Check input parameters
if ( ( i_response < 0 ) .or. ( i_response > 2 ) ) then
......
......@@ -8,6 +8,7 @@ module sca
LDA_we,LDA_ee,LDA_ww,LDA_rw,LDA_sww, LDA_s_ww_inv
integer :: LDA_ey,LDA_ye,LDA_dee
integer :: DESCA(9),DESCB(9),DESCZ(9),DESCC(9)
integer :: DESC_ye(9)
integer :: INFO,INFO_A,INFO_B,INFO_Z,INFO_C
real :: ORFAC
integer :: lwork_cooficient
......
......@@ -17,6 +17,7 @@ if(rank==0) write(*,*) 'output begings'
!======================================
!call print_starwall_m_ey
!======================================
!call print_starwall_m_ey2
!======================================
!call print_starwall_m_ee
......@@ -27,7 +28,8 @@ if(rank==0) write(*,*) 'output begings'
!======================================
!======================================
call print_starwall_response
!call print_starwall_response
call print_starwall_response2
!======================================
!======================================
......
......@@ -34,3 +34,54 @@ real :: num
if(rank==0) close(60)
end subroutine print_starwall_m_ey
subroutine print_starwall_m_ey2
use mpi_v
use resistive
use solv
use sca
!-----------------------------------------------------------------------
implicit none
include "mpif.h"
integer :: i_loc,j_loc, IC,JC
real :: num
integer INDXL2G
EXTERNAL INDXL2G
!-----------------------------------------------------------------------
call SCA_GRID(nd_bez,n_w)
CALL DESCINIT(DESCA, nd_bez, n_w, NB, NB, 0, 0, CONTEXT, LDA_ey, INFO_A )
if(INFO_A .NE. 0) then
write(6,*) "Something is wrong in print_starwall_m_ye CALL DESCINIT DESCA, INFO_A=",INFO_A
stop
endif
! if(rank==0) then
open(60,file='starwall_m_ey2',form="formatted",iostat=ier)
! write(60,'(2i8)') nd_bez,n_w
! endif
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
write(60,'(2i8,1pe16.8,1i8)') IC,JC,a_ey_loc(i_loc,j_loc), rank
END DO
END DO
!if(rank==0) close(60)
close(60)
end subroutine print_starwall_m_ey2
......@@ -364,3 +364,434 @@ end if !if (format_type == 'formatted' ) then
close(60)
end subroutine print_starwall_response
subroutine print_starwall_response2
use icontr
use contr_su
use solv
use tri_w
use tri_p
use coil2d
use sca
use mpi_v
use time
use resistive
!-----------------------------------------------------------------------
implicit none
include "mpif.h"
integer :: i,j,k,i_loc,j_loc
real :: num=0.
real,dimension(:,:),allocatable :: a_ye_print(:,:),a_ey_print(:,:),a_ee_print(:,:), &
s_ww_print(:,:),s_ww_inv_print(:,:)
integer :: IC, JC
integer INDXL2G
EXTERNAL INDXL2G
!-----------------------------------------------------------------------
! (Determine positions of wall triangle nodes)
if(rank==0) then
allocate( xyzpot_w(npot_w,3) )
do i = 1, ntri_w
do j = 1, 3
xyzpot_w(jpot_w(i,j),:) = (/ xw(i,j), yw(i,j), zw(i,j) /)
end do
end do
endif
if (format_type == 'formatted' ) then
if(rank==0) then
open(60, file='starwall-response2.dat', form=format_type, status='replace', action='write')
130 format(a)
131 format('#@intparam ',a24,99i12)
132 format('#@array ',a24,i12,a24,99i12)
133 format(4es24.16)
134 format(8i12)
write(60,130) '#@content STARWALL VACUUM RESPONSE DATA FOR THE JOREK CODE'
write(60,131) 'file_version', 1
write(60,131) 'n_bnd', N_bnd
write(60,131) 'nd_bez', nd_bez
write(60,131) 'ncoil', ncoil
write(60,131) 'npot_w', npot_w
write(60,131) 'n_w', n_w
write(60,131) 'ntri_w', ntri_w
write(60,131) 'n_tor', n_tor_jorek
write(60,132) 'i_tor', 1, 'int', n_tor_jorek, 0
write(60,134) i_tor_jorek(1:n_tor_jorek)
write(60,132) 'yy', 1, 'float', n_w, 0
write(60,133) 1./gamma(:)
write(60,132) 'ye', 2, 'float', n_w, nd_bez
endif
allocate(a_ye_print(n_w,nd_bez), stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate a_ye_print MYPROC_NUM=",MYPNUM
STOP
END IF
a_ye_print=0.
call SCA_GRID(n_w,nd_bez)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
a_ye_print(IC,JC) = a_ye_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, a_ye_print, n_w*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(a_ye_print, a_ye_print, n_w*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) write(60,133) a_ye_print(:,:)
deallocate (a_ye_print, a_ye_loc)
!=========================================================================================
allocate(a_ey_print(nd_bez,n_w),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate a_ey_print MYPROC_NUM=",MYPNUM
STOP
END IF
a_ey_print=0.
call SCA_GRID(nd_bez,n_w)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
a_ey_print(IC,JC) = a_ey_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, a_ey_print, nd_bez*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(a_ey_print, a_ey_print, nd_bez*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60,132) 'ey', 2, 'float', nd_bez, n_w
write(60,133) a_ey_print(:,:)
endif
deallocate (a_ey_print, a_ey_loc)
!=========================================================================================
allocate(a_ee_print(nd_bez,nd_bez),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate a_ee_print MYPROC_NUM=",MYPNUM
STOP
END IF
a_ee_print=0.
call SCA_GRID(nd_bez,nd_bez)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
a_ee_print(IC,JC) = a_ee_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, a_ee_print, nd_bez*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(a_ee_print, a_ee_print, nd_bez*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60,132) 'ee', 2, 'float', nd_bez, nd_bez
write(60,133) a_ee_print(:,:)
endif
deallocate (a_ee_print, a_ee_loc)
!=========================================================================================
allocate(s_ww_print(n_w,n_w),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate s_ww_print MYPROC_NUM=",MYPNUM
STOP
END IF
s_ww_print=0.
call SCA_GRID(n_w,n_w)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
s_ww_print(IC,JC) = s_ww_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, s_ww_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(s_ww_print, s_ww_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60,132) 's_ww', 2, 'float', n_w, n_w
write(60,133) s_ww_print(:,:)
endif
deallocate (s_ww_print, s_ww_loc)
!=========================================================================================
allocate(s_ww_inv_print(n_w,n_w),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate s_ww_inv_print MYPROC_NUM=",MYPNUM
STOP
END IF
s_ww_inv_print=0.
call SCA_GRID(n_w,n_w)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
s_ww_inv_print(IC,JC) = s_ww_inv_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, s_ww_inv_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(s_ww_inv_print, s_ww_inv_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60,132) 's_ww_inv', 2, 'float', n_w, n_w
write(60,133) s_ww_inv_print(:,:)
endif
deallocate (s_ww_inv_print, s_ww_inv_loc)
!=========================================================================================
if(rank==0) then
write(60,132) 'xyzpot_w', 2, 'float', npot_w, 3
write(60,133) xyzpot_w(:,:)
write(60,132) 'jpot_w', 2, 'int', ntri_w, 3
write(60,134) jpot_w(:,:)
deallocate(xyzpot_w)
endif
deallocate(jpot_w)
!=========================================================================================
else !if (format_type == 'unformatted' ) then
if(rank==0) then
open(60, file='starwall-response2.dat', form=format_type, status='replace', action='write')
char512='#@content STARWALL VACUUM RESPONSE DATA FOR THE JOREK CODE'
write(60) 42, 42.d0 !### for an elementary check in JOREK
write(60) char512
write(60) '#@intparam ', 'file_version ', 1
write(60) '#@intparam ', 'n_bnd ', N_bnd
write(60) '#@intparam ', 'nd_bez ', nd_bez
write(60) '#@intparam ', 'ncoil ', ncoil
write(60) '#@intparam ', 'npot_w ', npot_w
write(60) '#@intparam ', 'n_w ', n_w
write(60) '#@intparam ', 'ntri_w ', ntri_w
write(60) '#@intparam ', 'n_tor ', n_tor_jorek
write(60) '#@array ', 'i_tor ', 1, 'int ', n_tor_jorek, 0
write(60) i_tor_jorek(1:n_tor_jorek)
write(60) '#@array ', 'yy ', 1, 'float ', n_w, 0
write(60) 1.d0/gamma(:)
write(60) '#@array ', 'ye ', 2, 'float ', n_w, nd_bez
endif
allocate(a_ye_print(n_w,nd_bez), stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate a_ye_print MYPROC_NUM=",MYPNUM
STOP
END IF
a_ye_print=0.
call SCA_GRID(n_w,nd_bez)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
a_ye_print(IC,JC) = a_ye_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, a_ye_print, n_w*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(a_ye_print, a_ye_print, n_w*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) write(60) a_ye_print(:,:)
deallocate (a_ye_print, a_ye_loc)
!=========================================================================================
allocate(a_ey_print(nd_bez,n_w),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate a_ey_print MYPROC_NUM=",MYPNUM
STOP
END IF
a_ey_print=0.
call SCA_GRID(nd_bez,n_w)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
a_ey_print(IC,JC) = a_ey_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, a_ey_print, nd_bez*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(a_ey_print, a_ey_print, nd_bez*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60) '#@array ', 'ey ', 2, 'float ', nd_bez, n_w
write(60) a_ey_print(:,:)
endif
deallocate (a_ey_print, a_ey_loc)
!=========================================================================================
allocate(a_ee_print(nd_bez,nd_bez),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate a_ee_print MYPROC_NUM=",MYPNUM
STOP
END IF
a_ee_print=0.
call SCA_GRID(nd_bez,nd_bez)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
a_ee_print(IC,JC) = a_ee_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, a_ee_print, nd_bez*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(a_ee_print, a_ee_print, nd_bez*nd_bez, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60) '#@array ', 'ee ', 2, 'float ', nd_bez, nd_bez
write(60) a_ee_print(:,:)
endif
deallocate (a_ee_print, a_ee_loc)
!=========================================================================================
allocate(s_ww_print(n_w,n_w),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate s_ww_print MYPROC_NUM=",MYPNUM
STOP
END IF
s_ww_print=0.
call SCA_GRID(n_w,n_w)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
s_ww_print(IC,JC) = s_ww_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, s_ww_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(s_ww_print, s_ww_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60) '#@array ', 's_ww ', 2, 'float ', n_w, n_w
write(60) s_ww_print(:,:)
endif
deallocate (s_ww_print, s_ww_loc)
!=========================================================================================
allocate(s_ww_inv_print(n_w,n_w),stat=ier)
IF (IER /= 0) THEN
WRITE (*,*) "output.f90, can not allocate s_ww_inv_print MYPROC_NUM=",MYPNUM
STOP
END IF
s_ww_inv_print=0.
call SCA_GRID(n_w,n_w)
DO i_loc = 1,MP_A
IC= INDXL2G( i_loc, NB, MYROW, 0, NPROW)
DO j_loc = 1,NQ_A
JC= INDXL2G( j_loc, NB, MYCOL, 0, NPCOL)
s_ww_inv_print(IC,JC) = s_ww_inv_loc(i_loc,j_loc)
END DO
END DO
if(rank==0) then
call MPI_REDUCE(MPI_IN_PLACE, s_ww_inv_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
else
call MPI_REDUCE(s_ww_inv_print, s_ww_inv_print, n_w*n_w, &
MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, IER)
endif
if(rank==0) then
write(60) '#@array ', 's_ww_inv ', 2, 'float ', n_w, n_w
write(60) s_ww_inv_print(:,:)
endif
deallocate (s_ww_inv_print, s_ww_inv_loc)
!=========================================================================================
if(rank==0) then
write(60) '#@array ', 'xyzpot_w ', 2, 'float ', npot_w, 3
write(60) xyzpot_w(:,:)
write(60) '#@array ', 'jpot_w ', 2, 'int ', ntri_w, 3
write(60) jpot_w(:,:)
deallocate (xyzpot_w)
endif
deallocate (jpot_w)
!=========================================================================================
end if !if (format_type == 'formatted' ) then
close(60)
end subroutine print_starwall_response2
......@@ -80,6 +80,13 @@ if(rank==0) write(outp,*) ' n_w m_w rc_w rs_w zc_w
end do
if(rank==0) write(outp,*)
! --- Check input parameters
if ( ( nwu*nwv > 19600 ) ) then
write(outp,*) 'ERROR: npot_w (nwu*nwv) must have a value smaller than 19600 for current output version'
stop 1
end if
! Remarks:
! large Phi: independent variables
! small phi: current potential at each triangle vertex
......