gray/jetto interface reverted to F77 (removed use of modules), adding an intermediate F90 subroutine handling the reading of parameters/beams files. pec_init now uses optional input rho grid

This commit is contained in:
Lorenzo Figini 2015-11-24 16:36:20 +00:00
parent 46e36a5792
commit 59617c7a06
5 changed files with 242 additions and 57 deletions

136
Makefile.jetto Normal file
View File

@ -0,0 +1,136 @@
# Makefile for JETTO fortran 90 library 'gray'
# Author: Lorenzo Figini (figini@ifp.cnr.it)
#
# Derived from
# Makefile for JETTO fortran 77 library 'frantic'
# Author: Derek Harting (d.harting@fz-juelich.de)
# and
# Makefile for ITCequ Fortran 90 library - definitions
# P. Strand, elfps@chalmers.se
# G. Corrigan, gcor@jet.uk
# D. Harting, d.harting@fz-juelich.de
# Set the environment from the top-level Makefile of JETTO
# --------------------------------------------------------
#include ../include.mk
# Alternative search paths
vpath %.f90 src
JLIBDIR = ./
F90 = gfortran
DIRECTIVES = -DREVISION="'$(shell svnversion src)'"
# library name
# ------------
LIBNAME=$(JLIBDIR)/libgray.a
# List source and object files
# ----------------------------
FSRC=$(wildcard src/*.f90)
FOBJ0=$(subst src/,,$(FSRC))
FOBJ=$(FOBJ0:.f90=.o)
########################################
# Set up compiler specific environment #
########################################
# set default compiler
# --------------------
FC=$(F90)
# Set compiler debug and optimization flags
# -----------------------------------------
ifeq ("$(DBG)","")
# FFLAGS= $(AUTO) -O3 -Mpreprocess
FFLAGS=-O3 -cpp
else
# FFLAGS= $(AUTO) -O0 -Mpreprocess -g -Minform=inform -Mbounds -Mchkptr
FFLAGS=-Wall -g -finit-real=nan -ffpe-trap=invalid -fcheck=all -fbounds-check
endif
# Set include directories
# -----------------------
#INC=-I$(JCOMMON_STD)
INC=
#####################################
# Set up RULES specific environment #
#####################################
# suffixes
# --------
.SUFFIXES: .f90 .o .mod
# PHONY targets
# -------------
.PHONY: all clean realclean $(MASTER_RULES) gray
# default target
# --------------
all:
@($(MAKE) -C ../ all) || exit $$?
# Set rules passed to top-level Makefile
# --------------------------------------
$(MASTER_RULES):
@($(MAKE) -C ../ $@) || exit $$?
# rule for libgray.a
# --------------------
$(LIBNAME): $(FOBJ)
ar vr $(LIBNAME) $?
# synonym for libgray.a
# -----------------------
gray: $(LIBNAME)
# create rule for compiling f90 files
# -----------------------------------
%.o: %.f90
$(FC) -c $(DIRECTIVES) $(FFLAGS) $(INC) $<
# make 'realclean' option
# -----------------------
realclean:
rm -f *.o *.mod $(LIBNAME)
# make 'clean' option
# -------------------
clean:
rm -f *.o
# Dependencies
# ------------
main.o: const_and_precisions.o beams.o coreprofiles.o equilibrium.o \
graycore.o gray_params.o reflections.o
gray_jetto1beam.o: const_and_precisions.o beams.o graycore.o gray_params.o \
units.o
graycore.o: const_and_precisions.o beamdata.o beams.o coreprofiles.o \
dispersion.o eccd.o equilibrium.o errcodes.o gray_params.o \
pec.o polarization.o limiter.o units.o utils.o
beams.o: const_and_precisions.o dierckx.o reflections.o simplespline.o utils.o
beamdata.o: const_and_precisions.o gray_params.o
conical.o: const_and_precisions.o
coreprofiles.o: const_and_precisions.o dierckx.o gray_params.o simplespline.o \
utils.o
dierckx.o: const_and_precisions.o
dispersion.o: const_and_precisions.o eierf.o errcodes.o math.o quadpack.o
eccd.o: const_and_precisions.o conical.o dierckx.o errcodes.o magsurf_data.o \
numint.o
eierf.o: const_and_precisions.o
errcodes.o: const_and_precisions.o
gray_params.o: const_and_precisions.o utils.o
equilibrium.o: const_and_precisions.o dierckx.o limiter.o minpack.o \
reflections.o simplespline.o utils.o gray_params.o
magsurf_data.o: const_and_precisions.o gray_params.o equilibrium.o dierckx.o \
reflections.o simplespline.o units.o utils.o
math.o: const_and_precisions.o
minpack.o: const_and_precisions.o
numint.o: const_and_precisions.o
pec.o: const_and_precisions.o beamdata.o equilibrium.o gray_params.o \
magsurf_data.o utils.o
polarization.o: const_and_precisions.o
quadpack.o: const_and_precisions.o
reflections.o: const_and_precisions.o limiter.o utils.o
simplespline.o: const_and_precisions.o
utils.o: const_and_precisions.o

86
src/gray_jetto1beam.f90 Normal file
View File

@ -0,0 +1,86 @@
subroutine gray_jetto1beam(ijetto, mr, mz, r, z, psin, psia, rax, zax, &
nbnd, rbnd, zbnd, nrho, psrad, fpol, te, dne, zeff, qpsi, ibeam, &
powin, alphain, betain, dpdv, jcd, pabs, icd, ierr)
use const_and_precisions, only : wp_
use units, only : ucenr,usumm,uprj0,uwbm,udisp,ubres,ucnt,uoutr,ueq,uprfin, &
uflx,upec,uprm,ubeam
use gray_params, only : read_params,rtrparam_type,hcdparam_type,antctrl_type,&
eqparam_type,prfparam_type,outparam_type
use beams, only : read_beam2
use graycore, only : gray_main
implicit none
! arguments
integer, intent(in) :: ijetto, mr, mz, nbnd, nrho, ibeam
real(wp_), dimension(mr), intent(in) :: r
real(wp_), dimension(mz), intent(in) :: z
real(wp_), dimension(mr,mz), intent(in) :: psin
real(wp_), intent(in) :: psia, rax, zax, powin, alphain, betain
real(wp_), dimension(nbnd), intent(in) :: rbnd, zbnd
real(wp_), dimension(nrho), intent(in) :: psrad, fpol, te, dne, zeff, qpsi
real(wp_), dimension(nrho), intent(out) :: dpdv, jcd
real(wp_), intent(out) :: pabs, icd
integer, intent(out) :: ierr
! local variables
real(wp_), dimension(nrho) :: psinr
integer :: iox0
real(wp_) :: r0m,rvac,alpha0,beta0,psipol0,chipol0,p0mw,rwallm
real(wp_) :: fghz,x0,y0,z0,w1,w2,ri1,ri2,phiw,phir
real(wp_), dimension(5) :: rlim,zlim
logical, save :: firstcall=.true.
type(rtrparam_type) :: rtrp
type(hcdparam_type) :: hcdp
type(antctrl_type) :: antp
type(eqparam_type) :: eqp
type(prfparam_type) :: prfp
type(outparam_type) :: outp
! if first call read parameters from external file
if (firstcall) then
call read_params('gray.data',rtrp,hcdp,antp,eqp,rwallm,prfp,outp,uprm)
rwallm=r(1)
antp%filenm='graybeam.data'
eqp%filenm='JETTO'
eqp%iequil=ijetto+1
prfp%filenm='JETTO'
outp%ipec=1
firstcall=.false.
end if
rvac=rax
psinr=psrad
alpha0=alphain
beta0=betain
call read_beam2(antp%filenm,ibeam,alpha0,beta0,fghz,antp%iox,x0,y0,z0, &
w1,w2,ri1,ri2,phiw,phir,ubeam)
psipol0=antp%psi
chipol0=antp%chi
iox0=antp%iox
p0mw=powin*1.e-6_wp_
! set simple limiter
r0m=sqrt(x0**2+y0**2)*0.01_wp_
call range2rect(rwallm,max(r0m,r(mr)),z(1),z(mz),rlim,zlim)
! call main subroutine for the ibeam-th beam
call gray_main(r,z,psin,psia,psinr,fpol,qpsi,rvac,rax,zax,rbnd,zbnd,eqp, &
psrad,te,dne,zeff,prfp, rlim,zlim, p0mw,fghz,alpha0,beta0,(/x0,y0,z0/), &
w1,w2,ri1,ri2,phiw,phir,iox0,psipol0,chipol0, dpdv,jcd,pabs,icd,outp, &
rtrp,hcdp,ierr,sqrt(psrad))
! close output (debug) files
close(ucenr)
close(usumm)
close(uprj0)
close(uprj0+1)
close(uwbm)
close(udisp)
close(ubres)
close(ucnt)
close(uoutr)
close(ueq)
close(uprfin)
close(uflx)
close(upec)
end subroutine gray_jetto1beam

View File

@ -7,7 +7,7 @@ contains
subroutine gray_main(rv,zv,psin,psia,psinr,fpol,qpsi,rvac,rax,zax,rbnd,zbnd, &
eqp,psrad,terad,derad,zfc,prfp, rlim,zlim, &
p0,fghz,alpha0,beta0,xv0,w1,w2,ri1,ri2,phiw,phir,iox0, &
psipol0,chipol0, dpdv,jcd,pabs,icd, outp,rtrp,hcdp,ierr)
psipol0,chipol0, dpdv,jcd,pabs,icd, outp,rtrp,hcdp,ierr, rhout)
use const_and_precisions, only : zero, one
use coreprofiles, only : set_prfan, set_prfspl, temp, fzeff
use dispersion, only : expinit
@ -33,10 +33,10 @@ contains
type(rtrparam_type), intent(in) :: rtrp
type(hcdparam_type), intent(in) :: hcdp
real(wp_), dimension(:), allocatable, intent(in) :: psrad, terad, derad, zfc
real(wp_), dimension(:), allocatable, intent(in) :: rv, zv, psinr, fpol, qpsi
real(wp_), dimension(:), allocatable, intent(in) :: rbnd, zbnd, rlim, zlim
real(wp_), dimension(:,:), allocatable, intent(in) :: psin
real(wp_), dimension(:), intent(in) :: psrad, terad, derad, zfc
real(wp_), dimension(:), intent(in) :: rv, zv, psinr, fpol, qpsi
real(wp_), dimension(:), intent(in) :: rbnd, zbnd, rlim, zlim
real(wp_), dimension(:,:), intent(in) :: psin
real(wp_), intent(in) :: psia, rvac, rax, zax
integer, intent(in) :: iox0
real(wp_), intent(in) :: p0, fghz, psipol0, chipol0
@ -45,6 +45,7 @@ contains
real(wp_), intent(out) :: pabs,icd
real(wp_), dimension(:), intent(out) :: dpdv,jcd
real(wp_), dimension(:), intent(in), optional :: rhout
integer, intent(out) :: ierr
! local variables
@ -236,7 +237,7 @@ contains
if(nray > 1) call print_projxyzt(st,yw,1)
! compute power and current density profiles for all rays
call pec_init(ipec) !,sqrt(psinr))
call pec_init(ipec,rhout)
nnd=size(rhop_tab)
allocate(jphi(nnd),pins(nnd),currins(nnd))
call spec(psjki,ppabs,ccci,iiv,pabs,icd,dpdv,jphi,jcd,pins,currins)
@ -264,7 +265,7 @@ contains
! call dealloc_beam(yw,ypw,xc,du1,gri,ggri,psjki,ppabs,ccci, &
! tau0,alphaabs0,dids0,ccci0,p0jk,ext,eyt,iiv)
! call dealloc_pec
! deallocate(jphi,pins,currins)
deallocate(jphi,pins,currins)
! ======= free memory END ======
end subroutine gray_main

View File

@ -1,13 +1,13 @@
module units
! STANDARD
integer, parameter :: uprm = 2, ubeam = 97, uprf = 98, ueq = 99
integer, parameter :: uprfin = 55, uflx = 56, ubres = 70, ucnt = 71
integer, parameter :: uprj0 = 8, uwbm = 12, ucenr = 4, uoutr = 33
integer, parameter :: udisp = 17, upec = 48, usumm = 7
! integer, parameter :: uprm = 2, ubeam = 97, uprf = 98, ueq = 99
! integer, parameter :: uprfin = 55, uflx = 56, ubres = 70, ucnt = 71
! integer, parameter :: uprj0 = 8, uwbm = 12, ucenr = 4, uoutr = 33
! integer, parameter :: udisp = 17, upec = 48, usumm = 7
! JINTRAC
! integer, parameter :: uprm =602, ubeam =603, uprf =644, ueq =644
! integer, parameter :: uprfin =645, uflx =646, ubres =630, ucnt =631
! integer, parameter :: uprj0 =608, uwbm =612, ucenr =604, uoutr =633
! integer, parameter :: udisp =617, upec =648, usumm =607, uhead =638
integer, parameter :: uprm =602, ubeam =603, uprf =644, ueq =644
integer, parameter :: uprfin =645, uflx =646, ubres =630, ucnt =631
integer, parameter :: uprj0 =608, uwbm =612, ucenr =604, uoutr =633
integer, parameter :: udisp =617, upec =648, usumm =607, uhead =638
end module units

View File

@ -2,11 +2,6 @@
subroutine gray(ijetto, mr, mz, mrd, r, z, psin, psiax, psibnd,
. rax, zax, nbnd, rbnd, zbnd, nrho, psijet, f, te, dne, zeff,
. qsf, nbeam, powin, alphin, betain, dpdv, jcd, pec, icd, ierr)
use units, only : ucenr,usumm,uprj0,uprj0+1,uwbm,udisp,ubres,
. ucnt,uoutr,ueq,uprfin,uflx,upec
use gray_params, only : read_params
use beams, only : read_beam2
use graycore, only : gray_main
implicit none
! input arguments
integer ijetto, mr, mz, nbnd, nrho, nbeam
@ -21,10 +16,6 @@
! gray_main output arguments
real*8 dpdvloop(nrho), jcdloop(nrho), pecloop, icdloop
integer ierr
! local variables
real*8 rlim(5),zlim(5)
logical firstcall=.true.
save firstcall
! === input arguments ==================================================
!
@ -86,17 +77,6 @@
!
! ======================================================================
! if first call read parameters from external file
if (firstcall) then
call read_params('gray.data',rtrp,hcdp,antp,eqp,rwallm,
. prfp,outp,uprm)
antp%filenm='graybeam.data'
eqp%filenm='JETTO'
eqp%iequil=ijetto+1
prfp%filenm='JETTO'
firstcall=.false.
end if
! set output variables to 0
do i=1,nrho
dpdv(i) = 0.d0
@ -111,26 +91,12 @@
! read j-th beam properties from file
! and adjust alpha/beta if out of the allowed range
alpha0=alphin(j)
beta0=betain(j)
p0mw=powin(j)*1.d-6
call read_beam2(antp%filenm,j,alpha0,beta0,fghz,antp%iox,
. x0,y0,z0,w1,w2,ri1,ri2,phiw,phir,ubeam)
psipol0=antp%psi
chipol0=antp%chi
iox0=antp%iox
! set simple limiter
r0m=sqrt(x0**2+y0**2)*0.01d0
call range2rect(rwallm,max(r0m,rv(mr)),zv(1),zv(mz),rlim,zlim)
! call main subroutine for the j-th beam
subroutine gray_main(r,z,psin(1:mr,:),psibnd-psiax,
. psijet,-f,-qsf,rax,rax,zax,rbnd,zbnd,eqp,
. psijet,te,dne,zeff,prfp,rlim,zlim,
. p0mw,fghz,alpha0,beta0,(/x0,y0,z0/),
. w1,w2,ri1,ri2,phiw,phir,iox0,psipol0,chipol0,
. dpdvloop,jcdloop,pecpool,icdloop,outp,rtrp,hcdp,ierr)
subroutine gray_jetto1beam(ijetto, mr, mz, r, z, psin(1:mr,:),
. psibnd-psiax, rax, zax, nbnd, rbnd, zbnd, nrho, psijet, -f,
. te, dne, zeff, -qsf, j, powin(j), alphin(j), betain(j),
. dpdvloop, jcdloop, pecloop, icdloop, ierr)
! add contribution of j-th beam to the total
! adapting output data to JETTO convention on toroidal angle
@ -144,9 +110,5 @@
! end of loop over beams with power>0
end do
! close output (debug) files
close(ucenr,usumm,uprj0,uprj0+1,uwbm,udisp,ubres,ucnt,uoutr,
. ueq,uprfin,uflx,upec)
return
end