module io use itg_data use convert implicit none ! ! List of routines: ! ! wread ! Read gryffin .res file ! wpunch ! Write gryffin .res file ! wreadnc ! Read gryffin .nc file ! wpunchnc ! Write gryffin .nc file ! bwread ! Read gryffin .res file (unformatted) ! bwpunch ! Write gryffin .res file (unformatted) public ! Input/output routines for gryffin codes real, private, allocatable, dimension(:,:) :: ri1 real, private, allocatable, dimension(:,:,:) :: ri2 real, private, allocatable, dimension(:,:,:,:) :: ri3 real, private, allocatable, dimension(:,:,:,:,:) :: ri4 contains subroutine bwpunch(iou,ntim,ntimr,tim,dt1) ! ! Write/read results to/from the *.res file. ! The order of the *.res file is essentially given by the ! read/writes in this subroutine, and must be compatible between ! output subroutine wpunch.f and the input subroutine wread.f ! ! (c) Copyright 1991 to 1995 by Michael A. Beer, William D. Dorland, ! and Gregory W. Hammett. ALL RIGHTS RESERVED. ! ! subroutine arguments integer iou ! on entry, fortran i/o unit number of the *.res file integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. ! output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. ! local variables: integer i,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 real xxx ! For backwards compatibility with old cosine/sine coding ! (where the amplitude of the sine part is the negative of the ! amplitude of the imaginary part), reverse the sign of the imaginary ! parts of some of the fields. ! defaults values for backwards compatibility: xxx=0.0 tim=time(md) dt1=dt(md) ! start reading from or writing to the *.res file: write(iou)ld,ldb,kd,nmin,nmax,nd,lin,kdpass write(iou)nstp,ne,nfreq,ntotal,md,nspecies,ikx,iperiod write(iou)nparmom,nperpmom,nemom,iodd,iflr,iphi00,ifilter,iexp,igradon write(iou)(mmin(n),n=1,nd) write(iou)(mmax(n),n=1,nd) write(iou)tim,shr,qsf,epsn,eps,epse,nueeff,alpha, & etai,nuii,rmu1,xxx,xxx,tiovte,vy,vz,rmime,etae write(iou)dt0,dt1,x0,y0,z0,xp n1=nspecies nt1=ntim+1 nt2=ntim+ne ntim=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal ntimr=ntimr+ntotal ! eventually, put the following arrays in a more logical order: ! For now, keep backward compability. Put new arrays at end: allocate(ri4(2,ld,md,nd,n1)) call c2r(density, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) write(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !density write(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) allocate(ri3(2,ld,md,nd)) call c2r(potential, ri3) ri3(2,:,:,:) = - ri3(2,:,:,:) write(iou)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) !potential write(iou)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) deallocate (ri3) call c2r(t_par, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) write(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !t_par write(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(u_par, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) write(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !u_par write(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(t_perp, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) write(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !t_perp write(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(q_perp, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) write(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !q_perp write(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(q_par, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) write(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !q_par write(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) deallocate (ri4) ! The above replace 3rd and 4th perp moments, not yet used toroidally: ! write(iou)((((rps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! write(iou)((((rpc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! write(iou)((((sps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! write(iou)((((spc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) if (epse > 0.0) then allocate (ri3(2,kdpass,md,nd)) call c2r(e_density, ri3) write(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) ! e_density write(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(e_p, ri3) write(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_p write(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(e_r, ri3) write(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_r write(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(e_t, ri3) write(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_r write(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(phi_ba, ri3) write(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !phi_ba write(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) deallocate (ri3) allocate (ri3(2,ld,md,nd)) call c2r(phi_bk, ri3) write(iou)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) !phi_bk write(iou)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) call c2r(e_denk, ri3) write(iou)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) !e_denk write(iou)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) deallocate (ri3) allocate (ri2(2,md,nd)) call c2r(phiav_old, ri2) write(iou)((ri2(1,m,n),m=1,md),n=1,nd) ! phiav_old write(iou)((ri2(2,m,n),m=1,md),n=1,nd) deallocate (ri2) endif write(iou)(timo(i),i=nt1,nt2) !time history plot array write(iou)((wpfx(i,n),i=nt1,nt2),n=1,n1) !ion thermal flux write(iou)((fluxi(i,n),i=nt1,nt2),n=1,n1) !ion particle flux write(iou)(fluxe(i),i=nt1,nt2) !electron particle flux write(iou)(qfluxe(i),i=nt1,nt2) !electron thermal flux write(iou)(gamx(i),i=nt1,nt2) !total growth rate? write(iou)(((mgamx(m,n,i),m=1,md),n=1,nd),i=nt1,nt2) !gamma vs. mode ! write(iou)(wrhx(i),i=nt1,nt2) ! write(iou)(wrhy(i),i=nt1,nt2) ! write(iou)(wdux(i),i=nt1,nt2) ! write(iou)(wu2x(i),i=nt1,nt2) ! write(iou)(wp2x(i),i=nt1,nt2) ! write(iou)(wv2x(i),i=nt1,nt2) ! write(iou)(wlux(i),i=nt1,nt2) ! write(iou)(wdpx(i),i=nt1,nt2) ! write(iou)(wdvx(i),i=nt1,nt2) write(iou)(wenx(i),i=nt1,nt2) ! total energy ! write(iou)(wdfx(i),i=nt1,nt2) write(iou)(pcerr(i),i=nt1,nt2) ! energy conservation write(iou)(((wtif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode write(iou)(((wkif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !kin.en. vs. mode write(iou)(((wpif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !p^2 vs mode write(iou)(((psp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !phi*phi vs mode write(iou)(((wenrk(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode write(iou)(((dketa(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !eta drive v mode write(iou)(((dkpar(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !par damp. v mode write(iou)(((dktor(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor drive v mode write(iou)(((dktdp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor damp. v mode write(iou)(((wkups(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !wupsi v mode write(iou)(((dktot(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tot. drive v mode write(iou)((cx(i,n),i=nt1,nt2),n=1,4*nd) !phi(x) not aved. in y write(iou)((cy(i,m),i=nt1,nt2),m=1,4*md) !phi(y) not aved. in x write(iou)((cz(i,n),i=nt1,nt2),n=1,ld) !phi(z) not aved. in x,y write(iou)((czc(i,n),i=nt1,nt2),n=1,ld) !phi(0)phi(z) aved. in x,y write(iou)((czcn(i,n),i=nt1,nt2),n=1,ld) !phi(z)^2 aved. in x,y write(iou)((czn(i,n),i=nt1,nt2),n=1,ld) !ne(0)ne(z) aved. in x,y write(iou)((cznn(i,n),i=nt1,nt2),n=1,ld) !ne(z)^2 aved. in x,y write(iou)(ct(i),i=nt1,nt2) !phi(0)phi(t) aved. in x,y write(iou) ninterv write(iou)((phi00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. potential write(iou)((den00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. density write(iou)((upar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. velocity write(iou)((tpar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_par write(iou)((tperp00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_perp allocate (ri1(2,nt1:nt2)) call c2r(phih0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(uparh0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(denh0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(tparh0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(tperph0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(qparh0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(qperph0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(veh0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) call c2r(uparph0(nt1:nt2), ri1) write(iou)(ri1(1,i),i=nt1,nt2) write(iou)(ri1(2,i),i=nt1,nt2) deallocate (ri1) write(iou)((drt(id,i),id=1,14),i=nt1,nt2) !drive terms write(iou)(phirms(i),i=nt1,nt2) ! rms phi ! write(iou)(((wvif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! write(iou)(((wplx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) write(iou)(wakx(i),i=nt1,nt2) ! mode width in x write(iou)(waky(i),i=nt1,nt2) ! mode width in y write(iou)(wxsp(i),i=nt1,nt2) ! mode spread write(iou)(utor(i),i=nt1,nt2) ! toroidal flow write(iou)(upol(i),i=nt1,nt2) ! poloidal flow write(iou)(((grtmx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! amplitude allocate (ri3(2, ntf1:ntf2, md, nd)) call c2r(utim(ntf1:ntf2,:,:), ri3(:,ntf1:ntf2,:,:)) ! Phi(t) write(iou)(((ri3(1,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) !utim write(iou)(((ri3(2,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) deallocate (ri3) write(iou)(timf(i),i=ntf1,ntf2) ! time base for freq(t) and phi(t) write(iou)avgflag write(iou)date write(iou)note ! ! Write multi-species stuff at the end for backwards compatibility: ! write(iou)(rmass(n),n=1,nspecies) write(iou)(charge(n),n=1,nspecies) write(iou)(n_I(n),n=1,nspecies) write(iou)(Ln(n),n=1,nspecies) write(iou)(eta(n),n=1,nspecies) write(iou)(eta_par(n),n=1,nspecies) write(iou)(tau(n),n=1,nspecies) write(iou)(((phisq(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) write(iou)((((fluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) write(iou)((((qfluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) write(iou)(((fluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) write(iou)(((qfluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! more new variables tacked on at the end for upward compabitility: write(iou) icrit,ntimf write(iou)(time(m),m=1,md) write(iou)(dt(m),m=1,md) write(iou)(((rkperp2(l,m,n),l=1,ld),m=1,md),n=1,nd) write(iou) ((rdiff(m,n),m=1,md),n=1,nd) write(iou) chi_int write(iou) malias, nalias end subroutine bwpunch subroutine bwread(iou, ntim, ntimr, tim, dt1, nsteps, nfile) ! ! Write/read results to/from the *.res file in binary format. ! The order of the *.res file is essentially given by the ! read/writes in this subroutine, and must be compatible between ! output subroutine wpunch.f and the input subroutine wread.f ! ! (c) Copyright 1991 to 1995 by Michael A. Beer, William D. Dorland, ! and Gregory W. Hammett. ALL RIGHTS RESERVED. ! ! subroutine arguments integer iou ! on entry, fortran i/o unit number of the *.res file integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. integer nsteps ! (# of timesteps/nfreq) * (# of runs) ! for runs that have been restarted. integer nfile ! number of restart files ! output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. ! local variables: integer i,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 real xxx ! For backwards compatibility with old cosine/sine coding ! (where the amplitude of the sine part is the negative of the ! amplitude of the imaginary part), reverse the sign of the imaginary ! parts of some of the fields: ! defaults values for backwards compatibility: xxx=0.0 md=1 tim=time(md) dt1=dt(md) ! start reading from or writing to the *.res file: read(iou)ld,ldb,kd,nmin,nmax,nd,lin,kdpass read(iou)nstp,ne,nfreq,ntotal,md,nspecies,ikx,iperiod ! ! If mmin is not allocated, 'allocate' has not been called. ! if(.not.allocated(mmin)) then call allocate nsteps = max(nsteps, ntotal) call allocate_nez(nfile*ne, nsteps) endif read(iou)nparmom,nperpmom,nemom,iodd,iflr,iphi00,ifilter,iexp,igradon read(iou)(mmin(n),n=1,nd) read(iou)(mmax(n),n=1,nd) read(iou)tim,shr,qsf,epsn,eps,epse,nueeff,alpha, & etai,nuii,rmu1,xxx,xxx,tiovte,vy,vz,rmime,etae read(iou)dt0,dt1,x0,y0,z0,xp n1=nspecies nt1=ntim+1 nt2=ntim+ne ntim=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal ntimr=ntimr+ntotal write(*,*) ntf1,ntf2,ntimr,ntotal ! eventually, put the following arrays in a more logical order: ! For now, keep backward compability. Put new arrays at end: allocate (ri4(2,ld,md,nd,n1)) read(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !density read(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = - ri4(2,:,:,:,:) call r2c(density, ri4) allocate (ri3(2,ld,md,nd)) read(iou)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) !potential read(iou)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) ri3(2,:,:,:) = - ri3(2,:,:,:) call r2c(potential, ri3) deallocate (ri3) read(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !t_par read(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = - ri4(2,:,:,:,:) call r2c(t_par, ri4) read(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !u_par read(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = - ri4(2,:,:,:,:) call r2c(u_par, ri4) read(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !t_perp read(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = - ri4(2,:,:,:,:) call r2c(t_perp, ri4) read(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !q_perp read(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = - ri4(2,:,:,:,:) call r2c(q_perp, ri4) read(iou)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !q_par read(iou)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = - ri4(2,:,:,:,:) call r2c(q_par, ri4) deallocate (ri4) ! The above replace 3rd and 4th perp moments, not yet used toroidally: ! read(iou,2)((((rps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! read(iou,2)((((rpc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! read(iou,2)((((sps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! read(iou,2)((((spc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) if (epse > 0.0) then allocate (ri3(2,kdpass, md, nd)) read(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_density read(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_density, ri3) read(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_p read(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_p, ri3) read(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_r read(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_r, ri3) read(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !e_t read(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_t, ri3) read(iou)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) !phi_ba read(iou)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(phi_ba, ri3) deallocate (ri3) allocate (ri3(2,ld,md,nd)) read(iou)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) !phi_bk read(iou)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) call r2c(phi_bk, ri3) read(iou)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) !e_denk read(iou)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) call r2c(e_denk, ri3) deallocate (ri3) allocate (ri2(2,md,nd)) read(iou)((ri2(1,m,n),m=1,md),n=1,nd) !phiav_old read(iou)((ri2(2,m,n),m=1,md),n=1,nd) call r2c(phiav_old, ri2) deallocate (ri2) endif read(iou)(timo(i),i=nt1,nt2) !time history plot array read(iou)((wpfx(i,n),i=nt1,nt2),n=1,n1) !ion thermal flux read(iou)((fluxi(i,n),i=nt1,nt2),n=1,n1) !ion particle flux read(iou)(fluxe(i),i=nt1,nt2) !electron particle flux read(iou)(qfluxe(i),i=nt1,nt2) !electron thermal flux read(iou)(gamx(i),i=nt1,nt2) !total growth rate? read(iou)(((mgamx(m,n,i),m=1,md),n=1,nd),i=nt1,nt2) !gamma vs. mode ! read(iou)(wrhx(i),i=nt1,nt2) ! read(iou)(wrhy(i),i=nt1,nt2) ! read(iou)(wdux(i),i=nt1,nt2) ! read(iou)(wu2x(i),i=nt1,nt2) ! read(iou)(wp2x(i),i=nt1,nt2) ! read(iou)(wv2x(i),i=nt1,nt2) ! read(iou)(wlux(i),i=nt1,nt2) ! read(iou)(wdpx(i),i=nt1,nt2) ! read(iou)(wdvx(i),i=nt1,nt2) read(iou)(wenx(i),i=nt1,nt2) ! total energy ! read(iou)(wdfx(i),i=nt1,nt2) read(iou)(pcerr(i),i=nt1,nt2) ! energy conservation read(iou)(((wtif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode read(iou)(((wkif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !kin.en. vs. mode read(iou)(((wpif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !p^2 vs mode read(iou)(((psp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !phi*phi vs mode read(iou)(((wenrk(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode read(iou)(((dketa(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !eta drive v mode read(iou)(((dkpar(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !par damp. v mode read(iou)(((dktor(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor drive v mode read(iou)(((dktdp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor damp. v mode read(iou)(((wkups(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !wupsi v mode read(iou)(((dktot(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tot. drive v mode read(iou)((cx(i,n),i=nt1,nt2),n=1,4*nd) !phi(x) not aved. in y read(iou)((cy(i,m),i=nt1,nt2),m=1,4*md) !phi(y) not aved. in x read(iou)((cz(i,n),i=nt1,nt2),n=1,ld) !phi(z) not aved. in x,y read(iou)((czc(i,n),i=nt1,nt2),n=1,ld) !phi(0)phi(z) aved. in x,y read(iou)((czcn(i,n),i=nt1,nt2),n=1,ld) !phi(z)^2 aved. in x,y read(iou)((czn(i,n),i=nt1,nt2),n=1,ld) !ne(0)ne(z) aved. in x,y read(iou)((cznn(i,n),i=nt1,nt2),n=1,ld) !ne(z)^2 aved. in x,y read(iou)(ct(i),i=nt1,nt2) !phi(0)phi(t) aved. in x,y read(iou) ninterv read(iou)((phi00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. potential read(iou)((den00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. density read(iou)((upar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. velocity read(iou)((tpar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_par read(iou)((tperp00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_perp allocate (ri1(2,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(phih0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(uparh0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(denh0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(tparh0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(tperph0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(qparh0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(qperph0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(veh0, ri1(:,nt1:nt2)) read(iou)(ri1(1,i),i=nt1,nt2) read(iou)(ri1(2,i),i=nt1,nt2) call r2c(uparph0, ri1(:,nt1:nt2)) deallocate (ri1) read(iou)((drt(id,i),id=1,14),i=nt1,nt2) !drive terms read(iou)(phirms(i),i=nt1,nt2) ! rms phi ! read(iou)(((wvif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! read(iou)(((wplx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) read(iou)(wakx(i),i=nt1,nt2) ! mode width in x read(iou)(waky(i),i=nt1,nt2) ! mode width in y read(iou)(wxsp(i),i=nt1,nt2) ! mode spread read(iou)(utor(i),i=nt1,nt2) ! toroidal flow read(iou)(upol(i),i=nt1,nt2) ! poloidal flow read(iou)(((grtmx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! amplitude allocate (ri3(2,ntf1:ntf2,md,nd)) !utim read(iou)(((ri3(1,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) ! Phi(t) read(iou)(((ri3(2,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) call r2c(utim, ri3(:,ntf1:ntf2,:,:)) deallocate (ri3) read(iou)(timf(i),i=ntf1,ntf2) ! time base for freq(t) and phi(t) read(iou)avgflag read(iou)date read(iou)note ! ! Write multi-species stuff at the end for backwards compatibility: ! read(iou)(rmass(n),n=1,nspecies) read(iou)(charge(n),n=1,nspecies) read(iou)(n_I(n),n=1,nspecies) read(iou)(Ln(n),n=1,nspecies) read(iou)(eta(n),n=1,nspecies) read(iou)(eta_par(n),n=1,nspecies) read(iou)(tau(n),n=1,nspecies) read(iou)(((phisq(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) read(iou)((((fluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) read(iou)((((qfluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) read(iou)(((fluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) read(iou)(((qfluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! more new variables tacked on at the end for upward compabitility: read(iou) icrit,ntimf read(iou)(time(m),m=1,md) read(iou)(dt(m),m=1,md) read(iou)(((rkperp2(l,m,n),l=1,ld),m=1,md),n=1,nd) read(iou) ((rdiff(m,n),m=1,md),n=1,nd) read(iou) chi_int read(iou) malias, nalias end subroutine bwread subroutine wreadnc(ncid, ntim, ntimr, tim, dt1, nsteps, nfile) ! ! Read results from the *.nc file, and write them into the ! the itg common variables ! ! first read the dimensions and integers, then ! call separate subroutine to read rest of file include 'netcdf.inc' ! subroutine arguments integer ncid ! netCDF id number, from earlier nc_open call integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. integer nsteps ! total number of timesteps. This must be specified ! for runs that have been restarted. integer nfile ! number of restart files ! output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. ! local variables: character filename*80 ! name of netcdf file character yorn*1,varname*30 integer iwrite,ri,i,nt,k,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 integer i2(2),i3(3),i4(4),i5(5) real xxx integer emflag ! netCDF variables ! dimensions: integer :: status, nvars, ld_dim, md_dim, nd_dim, nspec_dim, & kdpass_dim, time_dim, nd4_dim, md4_dim, id_dim, timef_dim, & char30_dim, char80_dim integer :: ri_dim ! extra dimension for complex vars, 1=real 2=im ! variable tags: integer ld_id, ldb_id, kd_id, nmin_id, nmax_id integer nd_id, lin_id, kdpass_id integer nstp_id, ne_id, nfreq_id, ntotal_id, md_id integer nspecies_id, ikx_id, iperiod_id integer nparmom_id, nperpmom_id, nemom_id integer iodd_id, iflr_id, iphi00_id, ifilter_id integer iexp_id, igradon_id, mmin_id, mmax_id integer malias_id, nalias_id integer tim_id, shr_id, qsf_id, epsn_id, eps_id, epse_id integer nueeff_id, alpha_id, etai_id, nuii_id, rmu1_id integer tiovte_id, vy_id, vz_id, rmime_id, etae_id integer dt0_id, dt1_id, x0_id, y0_id, z0_id, xp_id integer density_id, potential_id, u_par_id, t_par_id, t_perp_id integer q_perp_id, q_par_id integer e_density_id, e_p_id, e_r_id, e_t_id, phi_ba_id integer phi_bk_id, e_denk_id, phiav_old_id integer semi_imp_id, beta_e_id ! EM integer n_e_id, u_e_id, tpar_e_id, tperp_e_id, apar_id ! EM integer timo_id, wpfx_id, fluxi_id, fluxe_id, qfluxe_id integer gamx_id, mgamx_id, wenx_id, pcerr_id integer wtif_id, wkif_id, wpif_id, psp_id, wenrk_id integer dketa_id, dkpar_id, dktor_id, dktdp_id, wkups_id integer dktot_id, cx_id, cy_id, cz_id, czc_id, czcn_id integer czn_id, cznn_id, ct_id, ninterv_id integer phi00_id, den00_id, upar00_id, tpar00_id, tperp00_id integer drt_id, phirms_id, wakx_id, waky_id, wxsp_id integer utor_id, upol_id, grtmx_id, utim_id, timf_id integer phih0_id,uparh0_id,denh0_id,tparh0_id,tperph0_id integer qparh0_id,qperph0_id,veh0_id,uparph0_id integer avgflag_id, date_id, note_id integer rmass_id, charge_id, n_I_id, Ln_id, eta_id integer eta_par_id, tau_id, phisq_id, fluximn_id integer qfluximn_id, fluxemn_id, qfluxemn_id integer icrit_id, ntimf_id, time_id, dt_id, rkperp2_id PRINT *, 'netCDF version ',NF_INQ_LIBVERS() ! open netCDF file ! 5/98 file is now opened by calling routine, ncid is id # ! filename=runname(1:lrunname)//'.nc' ! ! 10 STATUS = nf_open(filename, 0, NCID) ! read only ! IF (STATUS /= NF_NOERR) then ! PRINT *, NF_STRERROR(STATUS) ! print *, ' ' ! print *, 'Enter filename (with extension):' ! read *, filename ! goto 10 ! endif ! first determine if data file was produced by the electromagnetic ! version of the code by searching for the variable beta_e emflag=0 status=nf_inq_nvars(ncid,nvars) do i=1,nvars status=nf_inq_varname(ncid,i,varname) if (varname(1:6) == 'beta_e') then emflag=1 print *, 'Reading *.nc file from electromagnetic version' print *, ' ' endif enddo if (emflag == 0) print *, 'Reading *.nc file from electrostatic version' ! get dimensions id's status=nf_inq_dimid(ncid,'ld',ld_dim) status=nf_inq_dimid(ncid,'md',md_dim) status=nf_inq_dimid(ncid,'nd',nd_dim) status=nf_inq_dimid(ncid,'nspecies',nspec_dim) status=nf_inq_dimid(ncid,'ri',ri_dim) status=nf_inq_dimid(ncid,'time',time_dim) status=nf_inq_dimid(ncid,'timef',timef_dim) ! fast timescale status=nf_inq_dimid(ncid,'nd4',nd4_dim) status=nf_inq_dimid(ncid,'md4',md4_dim) status=nf_inq_dimid(ncid,'id',id_dim) ! drive term dim status=nf_inq_dimid(ncid,'char30',char30_dim) status=nf_inq_dimid(ncid,'char80',char80_dim) ! get variable id's status=nf_inq_varid(ncid,'ld',ld_id) status=nf_inq_varid(ncid,'ldb',ldb_id) status=nf_inq_varid(ncid,'kd',kd_id) status=nf_inq_varid(ncid,'nmin',nmin_id) status=nf_inq_varid(ncid,'nmax',nmax_id) status=nf_inq_varid(ncid,'nd',nd_id) status=nf_inq_varid(ncid,'lin',lin_id) status=nf_inq_varid(ncid,'kdpass',kdpass_id) status=nf_inq_varid(ncid,'nstp',nstp_id) status=nf_inq_varid(ncid,'ne',ne_id) status=nf_inq_varid(ncid,'nfreq',nfreq_id) status=nf_inq_varid(ncid,'ntotal',ntotal_id) status=nf_inq_varid(ncid,'md',md_id) status=nf_inq_varid(ncid,'nspecies',nspecies_id) status=nf_inq_varid(ncid,'ikx',ikx_id) status=nf_inq_varid(ncid,'iperiod',iperiod_id) status=nf_inq_varid(ncid,'nparmom',nparmom_id) status=nf_inq_varid(ncid,'nperpmom',nperpmom_id) status=nf_inq_varid(ncid,'nemom',nemom_id) status=nf_inq_varid(ncid,'iodd',iodd_id) status=nf_inq_varid(ncid,'iflr',iflr_id) status=nf_inq_varid(ncid,'iphi00',iphi00_id) status=nf_inq_varid(ncid,'ifilter',ifilter_id) status=nf_inq_varid(ncid,'iexp',iexp_id) status=nf_inq_varid(ncid,'igradon',igradon_id) status=nf_inq_varid(ncid,'malias',malias_id) status=nf_inq_varid(ncid,'nalias',nalias_id) status=nf_inq_varid(ncid,'mmin',mmin_id) status=nf_inq_varid(ncid,'mmax',mmax_id) status=nf_inq_varid(ncid,'tim',tim_id) status=nf_inq_varid(ncid,'shr',shr_id) status=nf_inq_varid(ncid,'qsf',qsf_id) status=nf_inq_varid(ncid,'epsn',epsn_id) status=nf_inq_varid(ncid,'eps',eps_id) status=nf_inq_varid(ncid,'epse',epse_id) status=nf_inq_varid(ncid,'nueeff',nueeff_id) status=nf_inq_varid(ncid,'alpha',alpha_id) status=nf_inq_varid(ncid,'etai',etai_id) status=nf_inq_varid(ncid,'nuii',nuii_id) status=nf_inq_varid(ncid,'rmu1',rmu1_id) status=nf_inq_varid(ncid,'tiovte',tiovte_id) status=nf_inq_varid(ncid,'vy',vy_id) status=nf_inq_varid(ncid,'vz',vz_id) status=nf_inq_varid(ncid,'rmime',rmime_id) status=nf_inq_varid(ncid,'etae',etae_id) status=nf_inq_varid(ncid,'dt0',dt0_id) status=nf_inq_varid(ncid,'dt1',dt1_id) status=nf_inq_varid(ncid,'x0',x0_id) status=nf_inq_varid(ncid,'y0',y0_id) status=nf_inq_varid(ncid,'z0',z0_id) status=nf_inq_varid(ncid,'xp',xp_id) status=nf_inq_varid(ncid,'potential',potential_id) status=nf_inq_varid(ncid,'density',density_id) status=nf_inq_varid(ncid,'u_par',u_par_id) status=nf_inq_varid(ncid,'t_par',t_par_id) status=nf_inq_varid(ncid,'t_perp',t_perp_id) status=nf_inq_varid(ncid,'q_par',q_par_id) status=nf_inq_varid(ncid,'q_perp',q_perp_id) status=nf_inq_varid(ncid,'phiav_old',phiav_old_id) status=nf_inq_varid(ncid,'timo',timo_id) status=nf_inq_varid(ncid,'wpfx',wpfx_id) status=nf_inq_varid(ncid,'fluxi',fluxi_id) status=nf_inq_varid(ncid,'fluxe',fluxe_id) status=nf_inq_varid(ncid,'qfluxe',qfluxe_id) status=nf_inq_varid(ncid,'gamx',gamx_id) status=nf_inq_varid(ncid,'mgamx',mgamx_id) status=nf_inq_varid(ncid,'wenx',wenx_id) status=nf_inq_varid(ncid,'pcerr',pcerr_id) status=nf_inq_varid(ncid,'wtif',wtif_id) status=nf_inq_varid(ncid,'wkif',wkif_id) status=nf_inq_varid(ncid,'wpif',wpif_id) status=nf_inq_varid(ncid,'psp',psp_id) status=nf_inq_varid(ncid,'wenrk',wenrk_id) status=nf_inq_varid(ncid,'dketa',dketa_id) status=nf_inq_varid(ncid,'dkpar',dkpar_id) status=nf_inq_varid(ncid,'dktor',dktor_id) status=nf_inq_varid(ncid,'dktdp',dktdp_id) status=nf_inq_varid(ncid,'wkups',wkups_id) status=nf_inq_varid(ncid,'dktot',dktot_id) status=nf_inq_varid(ncid,'cx',cx_id) status=nf_inq_varid(ncid,'cy',cy_id) status=nf_inq_varid(ncid,'cz',cz_id) status=nf_inq_varid(ncid,'czc',czc_id) status=nf_inq_varid(ncid,'czcn',czcn_id) status=nf_inq_varid(ncid,'czn',czn_id) status=nf_inq_varid(ncid,'cznn',cznn_id) status=nf_inq_varid(ncid,'ct',ct_id) status=nf_inq_varid(ncid,'ninterv',ninterv_id) status=nf_inq_varid(ncid,'phi00',phi00_id) status=nf_inq_varid(ncid,'den00',den00_id) status=nf_inq_varid(ncid,'upar00',upar00_id) status=nf_inq_varid(ncid,'tpar00',tpar00_id) status=nf_inq_varid(ncid,'tperp00',tperp00_id) status=nf_inq_varid(ncid,'phih0',phih0_id) status=nf_inq_varid(ncid,'uparh0',uparh0_id) status=nf_inq_varid(ncid,'denh0',denh0_id) status=nf_inq_varid(ncid,'tparh0',tparh0_id) status=nf_inq_varid(ncid,'tperph0',tperph0_id) status=nf_inq_varid(ncid,'qparh0',qparh0_id) status=nf_inq_varid(ncid,'qperph0',qperph0_id) status=nf_inq_varid(ncid,'veh0',veh0_id) status=nf_inq_varid(ncid,'uparph0',uparph0_id) status=nf_inq_varid(ncid,'drt',drt_id) status=nf_inq_varid(ncid,'phirms',phirms_id) status=nf_inq_varid(ncid,'wakx',wakx_id) status=nf_inq_varid(ncid,'waky',waky_id) status=nf_inq_varid(ncid,'wxsp',wxsp_id) status=nf_inq_varid(ncid,'utor',utor_id) status=nf_inq_varid(ncid,'upol',upol_id) status=nf_inq_varid(ncid,'grtmx',grtmx_id) status=nf_inq_varid(ncid,'utim',utim_id) status=nf_inq_varid(ncid,'timf',timf_id) status=nf_inq_varid(ncid,'avgflag',avgflag_id) status=nf_inq_varid(ncid,'date',date_id) status=nf_inq_varid(ncid,'note',note_id) status=nf_inq_varid(ncid,'rmass',rmass_id) status=nf_inq_varid(ncid,'charge',charge_id) status=nf_inq_varid(ncid,'n_I',n_I_id) status=nf_inq_varid(ncid,'Ln',Ln_id) status=nf_inq_varid(ncid,'eta',eta_id) status=nf_inq_varid(ncid,'eta_par',eta_par_id) status=nf_inq_varid(ncid,'tau',tau_id) status=nf_inq_varid(ncid,'phisq',phisq_id) status=nf_inq_varid(ncid,'fluximn',fluximn_id) status=nf_inq_varid(ncid,'qfluximn',qfluximn_id) status=nf_inq_varid(ncid,'fluxemn',fluxemn_id) status=nf_inq_varid(ncid,'qfluxemn',qfluxemn_id) status=nf_inq_varid(ncid,'icrit',icrit_id) status=nf_inq_varid(ncid,'ntimf',ntimf_id) status=nf_inq_varid(ncid,'time',time_id) status=nf_inq_varid(ncid,'dt',dt_id) status=nf_inq_varid(ncid,'rkperp2',rkperp2_id) ! read 0 and 1 dimensional variables status=nf_get_var_int(ncid,ld_id,ld) status=nf_get_var_int(ncid,ldb_id,ldb) status=nf_get_var_int(ncid,kd_id,kd) status=nf_get_var_int(ncid,nmin_id,nmin) status=nf_get_var_int(ncid,nmax_id,nmax) status=nf_get_var_int(ncid,nd_id,nd) status=nf_get_var_int(ncid,lin_id,lin) status=nf_get_var_int(ncid,kdpass_id,kdpass) status=nf_get_var_int(ncid,nstp_id,nstp) status=nf_get_var_int(ncid,ne_id,ne) status=nf_get_var_int(ncid,nfreq_id,nfreq) status=nf_get_var_int(ncid,ntotal_id,ntotal) status=nf_get_var_int(ncid,md_id,md) status=nf_get_var_int(ncid,nspecies_id,nspecies) status=nf_get_var_int(ncid,ikx_id,ikx) status=nf_get_var_int(ncid,iperiod_id,iperiod) status=nf_get_var_int(ncid,nparmom_id,nparmom) status=nf_get_var_int(ncid,nperpmom_id,nperpmom) status=nf_get_var_int(ncid,nemom_id,nemom) status=nf_get_var_int(ncid,iodd_id,iodd) status=nf_get_var_int(ncid,iflr_id,iflr) status=nf_get_var_int(ncid,iphi00_id,iphi00) status=nf_get_var_int(ncid,ifilter_id,ifilter) status=nf_get_var_int(ncid,iexp_id,iexp) status=nf_get_var_int(ncid,igradon_id,igradon) status=nf_get_var_int(ncid,malias_id,malias) status=nf_get_var_int(ncid,nalias_id,nalias) ! values for time indices: nt1=ntim+1 nt2=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal status=nf_get_var_double(ncid,tim_id,tim) status=nf_get_var_double(ncid,shr_id,shr) status=nf_get_var_double(ncid,qsf_id,qsf) status=nf_get_var_double(ncid,epsn_id,epsn) status=nf_get_var_double(ncid,eps_id,eps) status=nf_get_var_double(ncid,epse_id,epse) status=nf_get_var_double(ncid,nueeff_id,nueeff) status=nf_get_var_double(ncid,alpha_id,alpha) status=nf_get_var_double(ncid,etai_id,etai) status=nf_get_var_double(ncid,nuii_id,nuii) status=nf_get_var_double(ncid,rmu1_id,rmu1) status=nf_get_var_double(ncid,tiovte_id,tiovte) status=nf_get_var_double(ncid,vy_id,vy) status=nf_get_var_double(ncid,vz_id,vz) status=nf_get_var_double(ncid,rmime_id,rmime) status=nf_get_var_double(ncid,etae_id,etae) status=nf_get_var_double(ncid,dt0_id,dt0) status=nf_get_var_double(ncid,dt1_id,dt1) status=nf_get_var_double(ncid,x0_id,x0) status=nf_get_var_double(ncid,y0_id,y0) status=nf_get_var_double(ncid,z0_id,z0) status=nf_get_var_double(ncid,xp_id,xp) ! ! If mmin is not allocated, 'allocate' has not been called. ! if(.not.allocated(mmin)) then call allocate nsteps = max(nsteps, ntotal) ! ! This could be done more reliably by storing two new variables. ! One would count the total number of time steps taken to that point, ! and would replace nsteps, and the other would count the total number ! of times prstep was called, through all the continuation runs. This ! would replace nfile*ne here. ! call allocate_nez(nfile*ne, nsteps) endif status=nf_get_var_int(ncid,mmin_id,mmin) status=nf_get_var_int(ncid,mmax_id,mmax) status=nf_get_var_double(ncid,timo_id,timo(nt1:nt2)) status=nf_get_var_double(ncid,fluxe_id,fluxe(nt1:nt2)) status=nf_get_var_double(ncid,qfluxe_id,qfluxe(nt1:nt2)) status=nf_get_var_double(ncid,gamx_id,gamx(nt1:nt2)) status=nf_get_var_double(ncid,wenx_id,wenx(nt1:nt2)) status=nf_get_var_double(ncid,pcerr_id,pcerr(nt1:nt2)) status=nf_get_var_double(ncid,ct_id,ct(nt1:nt2)) status=nf_get_var_int(ncid,ninterv_id,ninterv) status=nf_get_var_double(ncid,phirms_id,phirms(nt1:nt2)) status=nf_get_var_double(ncid,wakx_id,wakx(nt1:nt2)) status=nf_get_var_double(ncid,waky_id,waky(nt1:nt2)) status=nf_get_var_double(ncid,wxsp_id,wxsp(nt1:nt2)) status=nf_get_var_double(ncid,utor_id,utor(nt1:nt2)) status=nf_get_var_double(ncid,upol_id,upol(nt1:nt2)) status=nf_get_var_double(ncid,timf_id,timf(ntf1:ntf2)) status=nf_get_var_double(ncid,avgflag_id,avgflag) status=nf_get_var_text(ncid,date_id,date) status=nf_get_var_text(ncid,note_id,note) status=nf_get_var_double(ncid,rmass_id,rmass) status=nf_get_var_double(ncid,charge_id,charge) status=nf_get_var_double(ncid,n_I_id,n_I) status=nf_get_var_double(ncid,Ln_id,Ln) status=nf_get_var_double(ncid,eta_id,eta) status=nf_get_var_double(ncid,eta_par_id,eta_par) status=nf_get_var_double(ncid,tau_id,tau) status=nf_get_var_int(ncid,icrit_id,icrit) status=nf_get_var_int(ncid,ntimf_id,ntimf) status=nf_get_var_double(ncid,time_id,time) status=nf_get_var_double(ncid,dt_id,dt) ! read multi-dimensional arrays n1=nspecies allocate (ri2(2,md,nd)) status=nf_get_var_double(ncid,phiav_old_id,ri2) call r2c (phiav_old, ri2) deallocate (ri2) allocate (ri3(2,ld,md,nd)) status=nf_get_var_double(ncid,potential_id,ri3) ri3(2,:,:,:) = -ri3(2,:,:,:) call r2c (potential, ri3) deallocate (ri3) status=nf_get_var_double(ncid,rkperp2_id,rkperp2(1:ld,1:md,1:nd)) allocate (ri4(2,ld,md,nd,n1)) status=nf_get_var_double(ncid,density_id,ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(density, ri4) status=nf_get_var_double(ncid,u_par_id,ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(u_par, ri4) status=nf_get_var_double(ncid,t_par_id,ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(t_par, ri4) status=nf_get_var_double(ncid,t_perp_id,ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(t_perp, ri4) status=nf_get_var_double(ncid,q_par_id,ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(q_par, ri4) status=nf_get_var_double(ncid,q_perp_id,ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(q_perp, ri4) deallocate (ri4) ! must be careful with the time variable to allow for calls with ! ntim or ntimr nonzero. the netCDF values will be 1...ne, but ! they must be read into nt1...nt2. status=nf_get_var_double(ncid,wpfx_id,wpfx(nt1:nt2,1:nspecies)) status=nf_get_var_double(ncid,fluxi_id,fluxi(nt1:nt2,1:nspecies)) status=nf_get_var_double(ncid,cx_id,cx(nt1:nt2,1:4*nd)) status=nf_get_var_double(ncid,cy_id,cy(nt1:nt2,1:4*md)) status=nf_get_var_double(ncid,cz_id,cz(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,czc_id,czc(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,czcn_id,czcn(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,czn_id,czn(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,cznn_id,cznn(nt1:nt2,1:ld)) status=nf_get_var_double(ncid,drt_id,drt(1:14,nt1:nt2)) allocate (ri2(2,nt1:nt2,nd)) status=nf_get_var_double(ncid,phi00_id,ri2) call r2c(phi00(nt1:nt2,:), ri2) status=nf_get_var_double(ncid,den00_id,ri2) call r2c(den00(nt1:nt2,:), ri2) status=nf_get_var_double(ncid,upar00_id,ri2) call r2c(upar00(nt1:nt2,:), ri2) status=nf_get_var_double(ncid,tpar00_id,ri2) call r2c(tpar00(nt1:nt2,:), ri2) status=nf_get_var_double(ncid,tperp00_id,ri2) call r2c(tperp00(nt1:nt2,:), ri2) deallocate (ri2) allocate (ri1(2,nt1:nt2)) status=nf_get_var_double(ncid,phih0_id,ri1) call r2c(phih0(nt1:nt2), ri1) status=nf_get_var_double(ncid,uparh0_id,ri1) call r2c(uparh0(nt1:nt2), ri1) status=nf_get_var_double(ncid,denh0_id,ri1) call r2c(denh0(nt1:nt2), ri1) status=nf_get_var_double(ncid,tparh0_id,ri1) call r2c(tparh0(nt1:nt2), ri1) status=nf_get_var_double(ncid,tperph0_id,ri1) call r2c(tperph0(nt1:nt2), ri1) status=nf_get_var_double(ncid,qparh0_id,ri1) call r2c(qparh0(nt1:nt2), ri1) status=nf_get_var_double(ncid,qperph0_id,ri1) call r2c(qperph0(nt1:nt2), ri1) status=nf_get_var_double(ncid,veh0_id,ri1) call r2c(veh0(nt1:nt2), ri1) status=nf_get_var_double(ncid,uparph0_id,ri1) call r2c(uparph0(nt1:nt2), ri1) deallocate (ri1) status=nf_get_var_double(ncid,mgamx_id,mgamx(1:md,1:nd,nt1:nt2)) status=nf_get_var_double(ncid,wtif_id,wtif(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wkif_id,wkif(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wpif_id,wpif(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,psp_id,psp(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wenrk_id,wenrk(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dketa_id,dketa(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dkpar_id,dkpar(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dktor_id,dktor(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dktdp_id,dktdp(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,wkups_id,wkups(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,dktot_id,dktot(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,grtmx_id,grtmx(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,phisq_id,phisq(nt1:nt2,1:md,1:nd)) if(epse > 0.) then status=nf_get_var_double(ncid,fluxemn_id,fluxemn(nt1:nt2,1:md,1:nd)) status=nf_get_var_double(ncid,qfluxemn_id,qfluxemn(nt1:nt2,1:md,1:nd)) endif status=nf_get_var_double(ncid,fluximn_id,fluximn(nt1:nt2,1:md,1:nd,1:n1)) status=nf_get_var_double(ncid,qfluximn_id,qfluximn(nt1:nt2,1:md,1:nd,1:n1)) allocate (ri3(2,ntotal,md,nd)) status=nf_get_var_double(ncid,utim_id,ri3) call r2c (utim(ntf1:ntf2,:,:), ri3) deallocate (ri3) ! set up bounce-avgd electrons only if they exist ! (avoids zero size dimension and saves disk space) if (epse > 0.0) then status=nf_inq_dimid(ncid,'kdpass',kdpass_dim) status=nf_inq_varid(ncid,'e_density',e_density_id) status=nf_inq_varid(ncid,'e_p',e_p_id) status=nf_inq_varid(ncid,'e_r',e_r_id) status=nf_inq_varid(ncid,'e_t',e_t_id) status=nf_inq_varid(ncid,'phi_ba',phi_ba_id) status=nf_inq_varid(ncid,'phi_bk',phi_bk_id) status=nf_inq_varid(ncid,'e_denk',e_denk_id) allocate (ri3(2, kdpass, md, nd)) status=nf_get_var_double(ncid,e_density_id,ri3) call r2c(e_density, ri3) status=nf_get_var_double(ncid,e_p_id,ri3) call r2c(e_p, ri3) status=nf_get_var_double(ncid,e_r_id,ri3) call r2c(e_r, ri3) status=nf_get_var_double(ncid,e_t_id,ri3) call r2c(e_t, ri3) status=nf_get_var_double(ncid,phi_ba_id,ri3) call r2c(phi_ba, ri3) deallocate (ri3) allocate (ri3(2,ld,md,nd)) status=nf_get_var_double(ncid,phi_bk_id,ri3) call r2c(phi_bk, ri3) status=nf_get_var_double(ncid,e_denk_id,ri3) call r2c(e_denk, ri3) deallocate (ri3) endif ! increment ntim and ntimr (for concatenation etc.) ntim=ntim+ne ntimr=ntimr+ntotal ! close file (now done by the calling routine) ! status = nf_close(ncid) ! if (status /= nf_noerr) stop end subroutine wreadnc !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc subroutine handle_err(status) include 'netcdf.inc' integer status if (status /= nf_noerr) then print *, nf_strerror(status) ! stop 'stopped' endif end subroutine handle_err subroutine wpunchnc(ncid,ntim,ntimr,tim,dt1) ! ! Write results to the *.nc file. ! include 'netcdf.inc' ! subroutine arguments integer ncid ! netCDF id number of the opened file integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. ! output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. ! local variables: character filename*80 ! name of netcdf file character yorn*1 integer i,k,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 real xxx integer ri ! netCDF variables integer status integer field_dim(4) ! dimensions of the complex field variables integer ion_dim(5) ! ion variables also have species # dimension integer e_dim(4) ! bounce-avgd electron dimensions integer av_dim(3) ! dimensions for flux-surface avgd complex integer flux_dim(2) ! time and nd integer mgamx_dim(3) ! md,nd,time integer wtif_dim(3) ! time,md,nd integer cx_dim(2),cy_dim(2),cz_dim(2),drt_dim(2) integer phi00_dim(3) ! ri,time,nd for flux-surface-avg time integer phih0_dim(2) ! ri,time integer utim_dim(4) integer fluximn_dim(4) integer lmn_dim(3) integer mn_dim(2) integer char30_dim, char80_dim integer ld_dim, md_dim, nd_dim, nspec_dim, kdpass_dim integer time_dim, nd4_dim, md4_dim, id_dim, timef_dim integer ri_dim ! extra dimension for complex vars, 1=real 2=im integer ld_id, ldb_id, kd_id, nmin_id, nmax_id ! variable tags integer nd_id, lin_id, kdpass_id integer nstp_id, ne_id, nfreq_id, ntotal_id, md_id integer nspecies_id, ikx_id, iperiod_id integer nparmom_id, nperpmom_id, nemom_id integer iodd_id, iflr_id, iphi00_id, ifilter_id integer iexp_id, igradon_id, mmin_id, mmax_id integer malias_id, nalias_id integer tim_id, shr_id, qsf_id, epsn_id, eps_id, epse_id integer nueeff_id, alpha_id, etai_id, nuii_id, rmu1_id integer tiovte_id, vy_id, vz_id, rmime_id, etae_id integer dt0_id, dt1_id, x0_id, y0_id, z0_id, xp_id integer density_id, potential_id, u_par_id, t_par_id, t_perp_id integer q_perp_id, q_par_id integer e_density_id, e_p_id, e_r_id, e_t_id, phi_ba_id integer phi_bk_id, e_denk_id, phiav_old_id integer timo_id, wpfx_id, fluxi_id, fluxe_id, qfluxe_id integer gamx_id, mgamx_id, wenx_id, pcerr_id integer wtif_id, wkif_id, wpif_id, psp_id, wenrk_id integer dketa_id, dkpar_id, dktor_id, dktdp_id, wkups_id integer dktot_id, cx_id, cy_id, cz_id, czc_id, czcn_id integer czn_id, cznn_id, ct_id, ninterv_id integer phi00_id, den00_id, upar00_id, tpar00_id, tperp00_id integer drt_id, phirms_id, wakx_id, waky_id, wxsp_id integer utor_id, upol_id, grtmx_id, utim_id, timf_id integer phih0_id,uparh0_id,denh0_id,tparh0_id,tperph0_id integer qparh0_id,qperph0_id,veh0_id,uparph0_id integer avgflag_id, date_id, note_id integer rmass_id, charge_id, n_I_id, Ln_id, eta_id integer eta_par_id, tau_id, phisq_id, fluximn_id integer qfluximn_id, fluxemn_id, qfluxemn_id integer icrit_id, ntimf_id, time_id, dt_id, rkperp2_id integer rdiff_id, chi_int_id ! small arrays for simple assign commands ! Is this variable not written?? It is used for simulations of particle noise. real rdiff_s(md,nd) ! defaults values for backwards compatibility: xxx=0.0 tim=time(md) dt1=dt(md) nt1=ntim+1 nt2=ntim+ne ntim=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal ntimr=ntimr+ntotal PRINT *, 'netCDF version ',NF_INQ_LIBVERS() ! create netCDF file ! 5/98 now done by calling routine ! filename=runname(1:lrunname)//'.nc' ! ! 10 STATUS = NF_CREATE(filename, 0, NCID) ! overwrites old ! IF (STATUS /= NF_NOERR) then ! PRINT *, NF_STRERROR(STATUS) ! print *, ' ' ! print *, 'retry? (y/n)' ! read *, yorn ! if (yorn == 'n') then ! stop ! else ! goto 10 ! endif ! endif ! create dimensions status=nf_def_dim(ncid,'ld',ld,ld_dim) status=nf_def_dim(ncid,'md',md,md_dim) status=nf_def_dim(ncid,'nd',nd,nd_dim) status=nf_def_dim(ncid,'nspecies',nspecies,nspec_dim) status=nf_def_dim(ncid,'ri',2,ri_dim) status=nf_def_dim(ncid,'time',ne,time_dim) status=nf_def_dim(ncid,'timef',ntotal,timef_dim) ! fast timescale status=nf_def_dim(ncid,'nd4',4*nd,nd4_dim) status=nf_def_dim(ncid,'md4',4*md,md4_dim) status=nf_def_dim(ncid,'id',14,id_dim) ! drive term dim status=nf_def_dim(ncid,'char30',30,char30_dim) status=nf_def_dim(ncid,'char80',80,char80_dim) ! set up bounce-avgd electrons only if they exist ! (avoids zero size dimension and saves disk space) if (epse > 0.0) then status=nf_def_dim(ncid,'kdpass',kdpass,kdpass_dim) endif field_dim(1)=ri_dim ! dimensions for phi, apar etc. field_dim(2)=ld_dim field_dim(3)=md_dim field_dim(4)=nd_dim ion_dim(1)=ri_dim ! dimensions for ion moments ion_dim(2)=ld_dim ion_dim(3)=md_dim ion_dim(4)=nd_dim ion_dim(5)=nspec_dim if (epse > 0.0) then e_dim(1)=ri_dim ! dimensions for bounce-avgd electron moments e_dim(2)=kdpass_dim e_dim(3)=md_dim e_dim(4)=nd_dim endif av_dim(1)=ri_dim ! dimensions for flux-surface avgd complex av_dim(2)=md_dim av_dim(3)=nd_dim flux_dim(1)=time_dim flux_dim(2)=nspec_dim mgamx_dim(1)=md_dim mgamx_dim(2)=nd_dim mgamx_dim(3)=time_dim wtif_dim(1)=time_dim wtif_dim(2)=md_dim wtif_dim(3)=nd_dim cx_dim(1)=time_dim cx_dim(2)=nd4_dim cy_dim(1)=time_dim cy_dim(2)=md4_dim cz_dim(1)=time_dim cz_dim(2)=ld_dim phi00_dim(1)=ri_dim phi00_dim(2)=time_dim phi00_dim(3)=nd_dim phih0_dim(1)=ri_dim phih0_dim(2)=time_dim drt_dim(1)=id_dim drt_dim(2)=time_dim utim_dim(1)=ri_dim utim_dim(2)=timef_dim utim_dim(3)=md_dim utim_dim(4)=nd_dim fluximn_dim(1)=time_dim fluximn_dim(2)=md_dim fluximn_dim(3)=nd_dim fluximn_dim(4)=nspec_dim lmn_dim(1)=ld_dim lmn_dim(2)=md_dim lmn_dim(3)=nd_dim mn_dim(1)=md_dim mn_dim(2)=nd_dim ! define variables status=nf_def_var(ncid,'ld',nf_int,0,0,ld_id) status=nf_def_var(ncid,'ldb',nf_int,0,0,ldb_id) status=nf_def_var(ncid,'kd',nf_int,0,0,kd_id) status=nf_def_var(ncid,'nmin',nf_int,0,0,nmin_id) status=nf_def_var(ncid,'nmax',nf_int,0,0,nmax_id) status=nf_def_var(ncid,'nd',nf_int,0,0,nd_id) status=nf_def_var(ncid,'lin',nf_int,0,0,lin_id) status=nf_def_var(ncid,'kdpass',nf_int,0,0,kdpass_id) status=nf_def_var(ncid,'nstp',nf_int,0,0,nstp_id) status=nf_def_var(ncid,'ne',nf_int,0,0,ne_id) status=nf_def_var(ncid,'nfreq',nf_int,0,0,nfreq_id) status=nf_def_var(ncid,'ntotal',nf_int,0,0,ntotal_id) status=nf_def_var(ncid,'md',nf_int,0,0,md_id) status=nf_def_var(ncid,'nspecies',nf_int,0,0,nspecies_id) status=nf_def_var(ncid,'ikx',nf_int,0,0,ikx_id) status=nf_def_var(ncid,'iperiod',nf_int,0,0,iperiod_id) status=nf_def_var(ncid,'nparmom',nf_int,0,0,nparmom_id) status=nf_def_var(ncid,'nperpmom',nf_int,0,0,nperpmom_id) status=nf_def_var(ncid,'nemom',nf_int,0,0,nemom_id) status=nf_def_var(ncid,'iodd',nf_int,0,0,iodd_id) status=nf_def_var(ncid,'iflr',nf_int,0,0,iflr_id) status=nf_def_var(ncid,'iphi00',nf_int,0,0,iphi00_id) status=nf_def_var(ncid,'ifilter',nf_int,0,0,ifilter_id) status=nf_def_var(ncid,'iexp',nf_int,0,0,iexp_id) status=nf_def_var(ncid,'igradon',nf_int,0,0,igradon_id) status=nf_def_var(ncid,'malias',nf_int,0,0,malias_id) status=nf_def_var(ncid,'nalias',nf_int,0,0,nalias_id) status=nf_def_var(ncid,'mmin',nf_int,1,nd_dim,mmin_id) status=nf_def_var(ncid,'mmax',nf_int,1,nd_dim,mmax_id) status=nf_def_var(ncid,'tim',nf_double,0,0,tim_id) status=nf_def_var(ncid,'shr',nf_double,0,0,shr_id) status=nf_def_var(ncid,'qsf',nf_double,0,0,qsf_id) status=nf_def_var(ncid,'epsn',nf_double,0,0,epsn_id) status=nf_def_var(ncid,'eps',nf_double,0,0,eps_id) status=nf_def_var(ncid,'epse',nf_double,0,0,epse_id) status=nf_def_var(ncid,'nueeff',nf_double,0,0,nueeff_id) status=nf_def_var(ncid,'alpha',nf_double,0,0,alpha_id) status=nf_def_var(ncid,'etai',nf_double,0,0,etai_id) status=nf_def_var(ncid,'nuii',nf_double,0,0,nuii_id) status=nf_def_var(ncid,'rmu1',nf_double,0,0,rmu1_id) status=nf_def_var(ncid,'tiovte',nf_double,0,0,tiovte_id) status=nf_def_var(ncid,'vy',nf_double,0,0,vy_id) status=nf_def_var(ncid,'vz',nf_double,0,0,vz_id) status=nf_def_var(ncid,'rmime',nf_double,0,0,rmime_id) status=nf_def_var(ncid,'etae',nf_double,0,0,etae_id) status=nf_def_var(ncid,'dt0',nf_double,0,0,dt0_id) status=nf_def_var(ncid,'dt1',nf_double,0,0,dt1_id) status=nf_def_var(ncid,'x0',nf_double,0,0,x0_id) status=nf_def_var(ncid,'y0',nf_double,0,0,y0_id) status=nf_def_var(ncid,'z0',nf_double,0,0,z0_id) status=nf_def_var(ncid,'xp',nf_double,0,0,xp_id) status=nf_def_var(ncid,'potential',nf_double,4,field_dim,potential_id) status=nf_def_var(ncid,'density',nf_double,5,ion_dim,density_id) status=nf_def_var(ncid,'u_par',nf_double,5,ion_dim,u_par_id) status=nf_def_var(ncid,'t_par',nf_double,5,ion_dim,t_par_id) status=nf_def_var(ncid,'t_perp',nf_double,5,ion_dim,t_perp_id) status=nf_def_var(ncid,'q_par',nf_double,5,ion_dim,q_par_id) status=nf_def_var(ncid,'q_perp',nf_double,5,ion_dim,q_perp_id) if (epse > 0.0) then status=nf_def_var(ncid,'e_density',nf_double,4,e_dim,e_density_id) status=nf_def_var(ncid,'e_p',nf_double,4,e_dim,e_p_id) status=nf_def_var(ncid,'e_r',nf_double,4,e_dim,e_r_id) status=nf_def_var(ncid,'e_t',nf_double,4,e_dim,e_t_id) status=nf_def_var(ncid,'phi_ba',nf_double,4,e_dim,phi_ba_id) status=nf_def_var(ncid,'phi_bk',nf_double,4,field_dim,phi_bk_id) status=nf_def_var(ncid,'e_denk',nf_double,4,field_dim,e_denk_id) endif status=nf_def_var(ncid,'phiav_old',nf_double,3,av_dim,phiav_old_id) status=nf_def_var(ncid,'timo',nf_double,1,time_dim,timo_id) status=nf_def_var(ncid,'wpfx',nf_double,2,flux_dim,wpfx_id) status=nf_def_var(ncid,'fluxi',nf_double,2,flux_dim,fluxi_id) status=nf_def_var(ncid,'fluxe',nf_double,1,time_dim,fluxe_id) status=nf_def_var(ncid,'qfluxe',nf_double,1,time_dim,qfluxe_id) status=nf_def_var(ncid,'gamx',nf_double,1,time_dim,gamx_id) status=nf_def_var(ncid,'mgamx',nf_double,3,mgamx_dim,mgamx_id) status=nf_def_var(ncid,'wenx',nf_double,1,time_dim,wenx_id) status=nf_def_var(ncid,'pcerr',nf_double,1,time_dim,pcerr_id) status=nf_def_var(ncid,'wtif',nf_double,3,wtif_dim,wtif_id) status=nf_def_var(ncid,'wkif',nf_double,3,wtif_dim,wkif_id) status=nf_def_var(ncid,'wpif',nf_double,3,wtif_dim,wpif_id) status=nf_def_var(ncid,'psp',nf_double,3,wtif_dim,psp_id) status=nf_def_var(ncid,'wenrk',nf_double,3,wtif_dim,wenrk_id) status=nf_def_var(ncid,'dketa',nf_double,3,wtif_dim,dketa_id) status=nf_def_var(ncid,'dkpar',nf_double,3,wtif_dim,dkpar_id) status=nf_def_var(ncid,'dktor',nf_double,3,wtif_dim,dktor_id) status=nf_def_var(ncid,'dktdp',nf_double,3,wtif_dim,dktdp_id) status=nf_def_var(ncid,'wkups',nf_double,3,wtif_dim,wkups_id) status=nf_def_var(ncid,'dktot',nf_double,3,wtif_dim,dktot_id) status=nf_def_var(ncid,'cx',nf_double,2,cx_dim,cx_id) status=nf_def_var(ncid,'cy',nf_double,2,cy_dim,cy_id) status=nf_def_var(ncid,'cz',nf_double,2,cz_dim,cz_id) status=nf_def_var(ncid,'czc',nf_double,2,cz_dim,czc_id) status=nf_def_var(ncid,'czcn',nf_double,2,cz_dim,czcn_id) status=nf_def_var(ncid,'czn',nf_double,2,cz_dim,czn_id) status=nf_def_var(ncid,'cznn',nf_double,2,cz_dim,cznn_id) status=nf_def_var(ncid,'ct',nf_double,1,time_dim,ct_id) status=nf_def_var(ncid,'ninterv',nf_int,0,0,ninterv_id) status=nf_def_var(ncid,'phi00',nf_double,3,phi00_dim,phi00_id) status=nf_def_var(ncid,'den00',nf_double,3,phi00_dim,den00_id) status=nf_def_var(ncid,'upar00',nf_double,3,phi00_dim,upar00_id) status=nf_def_var(ncid,'tpar00',nf_double,3,phi00_dim,tpar00_id) status=nf_def_var(ncid,'tperp00',nf_double,3,phi00_dim,tperp00_id) status=nf_def_var(ncid,'phih0',nf_double,2,phih0_dim,phih0_id) status=nf_def_var(ncid,'uparh0',nf_double,2,phih0_dim,uparh0_id) status=nf_def_var(ncid,'denh0',nf_double,2,phih0_dim,denh0_id) status=nf_def_var(ncid,'tparh0',nf_double,2,phih0_dim,tparh0_id) status=nf_def_var(ncid,'tperph0',nf_double,2,phih0_dim,tperph0_id) status=nf_def_var(ncid,'qparh0',nf_double,2,phih0_dim,qparh0_id) status=nf_def_var(ncid,'qperph0',nf_double,2,phih0_dim,qperph0_id) status=nf_def_var(ncid,'veh0',nf_double,2,phih0_dim,veh0_id) status=nf_def_var(ncid,'uparph0',nf_double,2,phih0_dim,uparph0_id) status=nf_def_var(ncid,'drt',nf_double,2,drt_dim,drt_id) status=nf_def_var(ncid,'phirms',nf_double,1,time_dim,phirms_id) status=nf_def_var(ncid,'wakx',nf_double,1,time_dim,wakx_id) status=nf_def_var(ncid,'waky',nf_double,1,time_dim,waky_id) status=nf_def_var(ncid,'wxsp',nf_double,1,time_dim,wxsp_id) status=nf_def_var(ncid,'utor',nf_double,1,time_dim,utor_id) status=nf_def_var(ncid,'upol',nf_double,1,time_dim,upol_id) status=nf_def_var(ncid,'grtmx',nf_double,3,wtif_dim,grtmx_id) status=nf_def_var(ncid,'utim',nf_double,4,utim_dim,utim_id) status=nf_def_var(ncid,'timf',nf_double,1,timef_dim,timf_id) status=nf_def_var(ncid,'avgflag',nf_double,0,0,avgflag_id) status=nf_def_var(ncid,'date',nf_char,1,char30_dim,date_id) status=nf_def_var(ncid,'note',nf_char,1,char80_dim,note_id) status=nf_def_var(ncid,'rmass',nf_double,1,nspec_dim,rmass_id) status=nf_def_var(ncid,'charge',nf_double,1,nspec_dim,charge_id) status=nf_def_var(ncid,'n_I',nf_double,1,nspec_dim,n_I_id) status=nf_def_var(ncid,'Ln',nf_double,1,nspec_dim,Ln_id) status=nf_def_var(ncid,'eta',nf_double,1,nspec_dim,eta_id) status=nf_def_var(ncid,'eta_par',nf_double,1,nspec_dim,eta_par_id) status=nf_def_var(ncid,'tau',nf_double,1,nspec_dim,tau_id) status=nf_def_var(ncid,'phisq',nf_double,3,wtif_dim,phisq_id) status=nf_def_var(ncid,'fluximn',nf_double,4,fluximn_dim,fluximn_id) status=nf_def_var(ncid,'qfluximn',nf_double,4,fluximn_dim,qfluximn_id) status=nf_def_var(ncid,'fluxemn',nf_double,3,wtif_dim,fluxemn_id) status=nf_def_var(ncid,'qfluxemn',nf_double,3,wtif_dim,qfluxemn_id) status=nf_def_var(ncid,'icrit',nf_int,0,0,icrit_id) status=nf_def_var(ncid,'ntimf',nf_int,0,0,ntimf_id) status=nf_def_var(ncid,'time',nf_double,1,md_dim,time_id) status=nf_def_var(ncid,'dt',nf_double,1,md_dim,dt_id) status=nf_def_var(ncid,'rkperp2',nf_double,3,lmn_dim,rkperp2_id) status=nf_def_var(ncid,'rdiff',nf_double,2,mn_dim,rdiff_id) status=nf_def_var(ncid,'chi_int',nf_double,0,0,chi_int,chi_int_id) ! add attributes status=NF_PUT_ATT_TEXT(NCID, NF_GLOBAL, 'title', 17, 'itgc restart data') status=nf_put_att_text(ncid,potential_id,'units',17,'T_i0/e rho_i/L_ne') status=nf_put_att_text(ncid,potential_id,'long_name',23,'Electrostatic Potential') status=nf_put_att_text(ncid,potential_id,'idl_name',5,'!7U!6') status=nf_put_att_text(ncid,density_id,'units',14,'n_0 rho_i/L_ne') status=nf_put_att_text(ncid,density_id,'long_name',11,'Ion Density') status=nf_put_att_text(ncid,density_id,'idl_name',8,'!6n!Di!N') status=nf_put_att_text(ncid,u_par_id,'units',15,'v_ti rho_i/L_ne') status=nf_put_att_text(ncid,u_par_id,'long_name',21,'Ion Parallel Velocity') status=nf_put_att_text(ncid,u_par_id,'idl_name',15,'!6u!D!9#!6!Di!N') status=nf_put_att_text(ncid,t_par_id,'units',15,'T_i0 rho_i/L_ne') status=nf_put_att_text(ncid,t_par_id,'long_name',24,'Ion Parallel Temperature') status=nf_put_att_text(ncid,t_par_id,'idl_name',15,'!6T!D!9#!6!Di!N') status=nf_put_att_text(ncid,t_perp_id,'units',15,'T_i0 rho_i/L_ne') status=nf_put_att_text(ncid,t_perp_id,'long_name',29,'Ion Perpendicular Temperature') status=nf_put_att_text(ncid,t_perp_id,'idl_name',15,'!6T!D!9x!6!Di!N') status=nf_put_att_text(ncid,q_par_id,'units',25,'n_0 m_i v_ti^3 rho_i/L_ne') status=nf_put_att_text(ncid,q_par_id,'long_name',22,'Ion Parallel Heat Flow') status=nf_put_att_text(ncid,q_par_id,'idl_name',15,'!6q!D!9#!6!Di!N') status=nf_put_att_text(ncid,q_perp_id,'units',25,'n_0 m_i v_ti^3 rho_i/L_ne') status=nf_put_att_text(ncid,q_perp_id,'long_name',27,'Ion Perpendicular Heat Flow') status=nf_put_att_text(ncid,q_perp_id,'idl_name',15,'!6q!D!9x!6!Di!N') status=nf_put_att_text(ncid,wpfx_id,'long_name',13,'Ion Heat Flux') status=nf_put_att_text(ncid,fluxi_id,'long_name',17,'Ion Particle Flux') status=nf_put_att_text(ncid,fluxe_id,'long_name',22,'Electron Particle Flux') status=nf_put_att_text(ncid,qfluxe_id,'long_name',18,'Electron Heat Flux') status=nf_put_att_text(ncid,gamx_id,'long_name',18,'Energy Growth Rate') status=nf_put_att_text(ncid,wenx_id,'long_name',12,'Total Energy') status=nf_put_att_text(ncid,ct_id,'long_name',29,'Time correlation of Potential') status=nf_put_att_text(ncid,phirms_id,'long_name',13,'RMS Potential') status=nf_put_att_text(ncid,wakx_id,'long_name',17,'Radial Mode Width') status=nf_put_att_text(ncid,waky_id,'long_name',19,'Poloidal Mode Width') status=nf_put_att_text(ncid,wxsp_id,'long_name',11,'Mode Spread') status=nf_put_att_text(ncid,utor_id,'long_name',13,'Toroidal Flow') status=nf_put_att_text(ncid,upol_id,'long_name',13,'Poloidal Flow') ! write variables status=nf_enddef(ncid) ! out of definition mode status=nf_put_var_int(ncid,ld_id,ld) status=nf_put_var_int(ncid,ldb_id,ldb) status=nf_put_var_int(ncid,kd_id,kd) status=nf_put_var_int(ncid,nmin_id,nmin) status=nf_put_var_int(ncid,nmax_id,nmax) status=nf_put_var_int(ncid,nd_id,nd) status=nf_put_var_int(ncid,lin_id,lin) status=nf_put_var_int(ncid,kdpass_id,kdpass) status=nf_put_var_int(ncid,nstp_id,nstp) status=nf_put_var_int(ncid,ne_id,ne) status=nf_put_var_int(ncid,nfreq_id,nfreq) status=nf_put_var_int(ncid,ntotal_id,ntotal) status=nf_put_var_int(ncid,md_id,md) status=nf_put_var_int(ncid,nspecies_id,nspecies) status=nf_put_var_int(ncid,ikx_id,ikx) status=nf_put_var_int(ncid,iperiod_id,iperiod) status=nf_put_var_int(ncid,nparmom_id,nparmom) status=nf_put_var_int(ncid,nperpmom_id,nperpmom) status=nf_put_var_int(ncid,nemom_id,nemom) status=nf_put_var_int(ncid,iodd_id,iodd) status=nf_put_var_int(ncid,iflr_id,iflr) status=nf_put_var_int(ncid,iphi00_id,iphi00) status=nf_put_var_int(ncid,ifilter_id,ifilter) status=nf_put_var_int(ncid,iexp_id,iexp) status=nf_put_var_int(ncid,igradon_id,igradon) status=nf_put_var_int(ncid,malias_id,malias) status=nf_put_var_int(ncid,nalias_id,nalias) status=nf_put_var_int(ncid,mmin_id,mmin) status=nf_put_var_int(ncid,mmax_id,mmax) status=nf_put_var_double(ncid,tim_id,tim) status=nf_put_var_double(ncid,shr_id,shr) status=nf_put_var_double(ncid,qsf_id,qsf) status=nf_put_var_double(ncid,epsn_id,epsn) status=nf_put_var_double(ncid,eps_id,eps) status=nf_put_var_double(ncid,epse_id,epse) status=nf_put_var_double(ncid,nueeff_id,nueeff) status=nf_put_var_double(ncid,alpha_id,alpha) status=nf_put_var_double(ncid,etai_id,etai) status=nf_put_var_double(ncid,nuii_id,nuii) status=nf_put_var_double(ncid,rmu1_id,rmu1) status=nf_put_var_double(ncid,tiovte_id,tiovte) status=nf_put_var_double(ncid,vy_id,vy) status=nf_put_var_double(ncid,vz_id,vz) status=nf_put_var_double(ncid,rmime_id,rmime) status=nf_put_var_double(ncid,etae_id,etae) status=nf_put_var_double(ncid,dt0_id,dt0) status=nf_put_var_double(ncid,dt1_id,dt1) status=nf_put_var_double(ncid,x0_id,x0) status=nf_put_var_double(ncid,y0_id,y0) status=nf_put_var_double(ncid,z0_id,z0) status=nf_put_var_double(ncid,xp_id,xp) allocate (ri3(2,ld,md,nd)) call c2r(potential, ri3) ri3(2,:,:,:) = - ri3(2,:,:,:) status=nf_put_var_double(ncid,potential_id,ri3) deallocate (ri3) allocate (ri4(2,ld,md,nd,nspecies)) call c2r(density, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) status=nf_put_var_double(ncid,density_id,ri4) call c2r(u_par, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) status=nf_put_var_double(ncid,u_par_id,ri4) call c2r(t_par, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) status=nf_put_var_double(ncid,t_par_id,ri4) call c2r(t_perp, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) status=nf_put_var_double(ncid,t_perp_id,ri4) call c2r(q_par, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) status=nf_put_var_double(ncid,q_par_id,ri4) call c2r(q_perp, ri4) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) status=nf_put_var_double(ncid,q_perp_id,ri4) deallocate (ri4) if (epse > 0.0) then allocate (ri3(2,kdpass, md, nd)) call c2r(e_density, ri3) status=nf_put_var_double(ncid,e_density_id,ri3) call c2r(e_p, ri3) status=nf_put_var_double(ncid,e_p_id,ri3) call c2r(e_r, ri3) status=nf_put_var_double(ncid,e_r_id,ri3) call c2r(e_t, ri3) status=nf_put_var_double(ncid,e_t_id,ri3) call c2r(phi_ba, ri3) status=nf_put_var_double(ncid,phi_ba_id,ri3) deallocate (ri3) allocate (ri3(2, ld, md, nd)) call c2r(phi_bk, ri3) status=nf_put_var_double(ncid,phi_bk_id,ri3) call c2r(e_denk, ri3) status=nf_put_var_double(ncid,e_denk_id,ri3) deallocate (ri3) endif allocate (ri2(2,md,nd)) call c2r(phiav_old, ri2) status=nf_put_var_double(ncid,phiav_old_id,ri2) deallocate (ri2) status=nf_put_var_double(ncid,timo_id,timo(nt1)) status=nf_put_var_double(ncid,wpfx_id,wpfx) status=nf_put_var_double(ncid,fluxi_id,fluxi) status=nf_put_var_double(ncid,fluxe_id,fluxe(nt1)) status=nf_put_var_double(ncid,qfluxe_id,qfluxe(nt1)) status=nf_put_var_double(ncid,gamx_id,gamx(nt1)) status=nf_put_var_double(ncid,mgamx_id,mgamx) status=nf_put_var_double(ncid,wenx_id,wenx(nt1)) status=nf_put_var_double(ncid,pcerr_id,pcerr(nt1)) status=nf_put_var_double(ncid,wtif_id,wtif) status=nf_put_var_double(ncid,wkif_id,wkif) status=nf_put_var_double(ncid,wpif_id,wpif) status=nf_put_var_double(ncid,psp_id,psp) status=nf_put_var_double(ncid,wenrk_id,wenrk) status=nf_put_var_double(ncid,dketa_id,dketa) status=nf_put_var_double(ncid,dkpar_id,dkpar) status=nf_put_var_double(ncid,dktor_id,dktor) status=nf_put_var_double(ncid,dktdp_id,dktdp) status=nf_put_var_double(ncid,wkups_id,wkups) status=nf_put_var_double(ncid,dktot_id,dktot) status=nf_put_var_double(ncid,cx_id,cx) status=nf_put_var_double(ncid,cy_id,cy) status=nf_put_var_double(ncid,cz_id,cz) status=nf_put_var_double(ncid,czc_id,czc) status=nf_put_var_double(ncid,czcn_id,czcn) status=nf_put_var_double(ncid,czn_id,czn) status=nf_put_var_double(ncid,cznn_id,cznn) status=nf_put_var_double(ncid,ct_id,ct(nt1)) status=nf_put_var_int(ncid,ninterv_id,ninterv) allocate (ri2(2,ne,nd)) call c2r(phi00(nt1:nt2,:), ri2) status=nf_put_var_double(ncid,phi00_id,ri2) call c2r(den00(nt1:nt2,:), ri2) status=nf_put_var_double(ncid,den00_id,ri2) call c2r(upar00(nt1:nt2,:), ri2) status=nf_put_var_double(ncid,upar00_id,ri2) call c2r(tpar00(nt1:nt2,:), ri2) status=nf_put_var_double(ncid,tpar00_id,ri2) call c2r(tperp00(nt1:nt2,:), ri2) status=nf_put_var_double(ncid,tperp00_id,ri2) deallocate (ri2) allocate (ri1(2,ne)) call c2r(phih0(nt1:nt2), ri1) status=nf_put_var_double(ncid,phih0_id,ri1) call c2r(uparh0(nt1:nt2), ri1) status=nf_put_var_double(ncid,uparh0_id,ri1) call c2r(denh0(nt1:nt2), ri1) status=nf_put_var_double(ncid,denh0_id,ri1) call c2r(tparh0(nt1:nt2), ri1) status=nf_put_var_double(ncid,tparh0_id,ri1) call c2r(tperph0(nt1:nt2), ri1) status=nf_put_var_double(ncid,tperph0_id,ri1) call c2r(qparh0(nt1:nt2), ri1) status=nf_put_var_double(ncid,qparh0_id,ri1) call c2r(qperph0(nt1:nt2), ri1) status=nf_put_var_double(ncid,qperph0_id,ri1) call c2r(veh0(nt1:nt2), ri1) status=nf_put_var_double(ncid,veh0_id,ri1) call c2r(uparph0(nt1:nt2), ri1) status=nf_put_var_double(ncid,uparph0_id,ri1) deallocate (ri1) status=nf_put_var_double(ncid,drt_id,drt) status=nf_put_var_double(ncid,phirms_id,phirms(nt1)) status=nf_put_var_double(ncid,wakx_id,wakx(nt1)) status=nf_put_var_double(ncid,waky_id,waky(nt1)) status=nf_put_var_double(ncid,wxsp_id,wxsp(nt1)) status=nf_put_var_double(ncid,utor_id,utor(nt1)) status=nf_put_var_double(ncid,upol_id,upol(nt1)) status=nf_put_var_double(ncid,grtmx_id,grtmx) allocate (ri3(2, ntotal, md, nd)) call c2r(utim(ntf1:ntf2,:,:), ri3) status=nf_put_var_double(ncid,utim_id,ri3) deallocate (ri3) status=nf_put_var_double(ncid,timf_id,timf(ntf1)) status=nf_put_var_double(ncid,avgflag_id,avgflag) status=nf_put_var_text(ncid,date_id,date) status=nf_put_var_text(ncid,note_id,note) status=nf_put_var_double(ncid,rmass_id,rmass) status=nf_put_var_double(ncid,charge_id,charge) status=nf_put_var_double(ncid,n_I_id,n_I) status=nf_put_var_double(ncid,Ln_id,Ln) status=nf_put_var_double(ncid,eta_id,eta) status=nf_put_var_double(ncid,eta_par_id,eta_par) status=nf_put_var_double(ncid,tau_id,tau) status=nf_put_var_double(ncid,phisq_id,phisq) status=nf_put_var_double(ncid,fluximn_id,fluximn) status=nf_put_var_double(ncid,qfluximn_id,qfluximn) if(epse > 0.) then status=nf_put_var_double(ncid,fluxemn_id,fluxemn) status=nf_put_var_double(ncid,qfluxemn_id,qfluxemn) endif status=nf_put_var_int(ncid,icrit_id,icrit) status=nf_put_var_int(ncid,ntimf_id,ntimf) status=nf_put_var_double(ncid,time_id,time) status=nf_put_var_double(ncid,dt_id,dt) status=nf_put_var_double(ncid,rkperp2_id,rkperp2) ! close file (now done by calling routine) ! status = nf_close(ncid) ! end subroutine wpunchnc subroutine wread(iou,ntim,ntimr,tim,dt1, nsteps, nfile) ! ! Write/read results to/from the *.res file. ! The format of the *.res file is essentially given by the ! read/writes in this subroutine, and must be compatible between ! output subroutine wpunch.f and the input subroutine wread.f ! ! (c) Copyright 1991 to 1995 by Michael A. Beer, William D. Dorland, ! and Gregory W. Hammett. ALL RIGHTS RESERVED. ! ! subroutine arguments integer iou ! on entry, fortran i/o unit number of the *.res file integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. integer nsteps ! (# of timesteps/nfreq) * (# of runs) ! for runs that have been restarted. integer nfile ! number of restart files ! output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. ! local variables: integer i,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 real xxx ! For backwards compatibility with old cosine/sine coding ! (where the amplitude of the sine part is the negative of the ! amplitude of the imaginary part), reverse the sign of the imaginary ! parts of some of the fields. ! defaults values for backwards compatibility: xxx=0.0 md=1 tim=time(md) dt1=dt(md) ! start reading from or writing to the *.res file: read(iou,1)ld,ldb,kd,nmin,nmax,nd,lin,kdpass read(iou,1)nstp,ne,nfreq,ntotal,md,nspecies,ikx,iperiod read(iou,1)nparmom,nperpmom,nemom,iodd,iflr,iphi00,ifilter,iexp,igradon ! ! If mmin is not allocated, 'allocate' has not been called. ! if(.not.allocated(mmin)) then call allocate nsteps = max(nsteps, ntotal) call allocate_nez(nfile*ne, nsteps) endif read(iou,1)(mmin(n),n=1,nd) read(iou,1)(mmax(n),n=1,nd) read(iou,2)tim,shr,qsf,epsn,eps,epse,nueeff,alpha, & etai,nuii,rmu1,xxx,xxx,tiovte,vy,vz,rmime,etae read(iou,2)dt0,dt1,x0,y0,z0,xp n1=nspecies nt1=ntim+1 nt2=ntim+ne ntim=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal ntimr=ntimr+ntotal ! eventually, put the following arrays in a more logical order: ! For now, keep backward compability. Put new arrays at end: allocate(ri4(2,ld,md,nd,n1)) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(density, ri4) allocate (ri3(2,ld,md,nd)) read(iou,2)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) read(iou,2)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) ri3(2,:,:,:) = -ri3(2,:,:,:) call r2c(potential, ri3) deallocate (ri3) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) !t_par read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(t_par, ri4) ! redundant write, just to keep other arrays in the same position ! in the file for backwards compatability, for now: read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(u_par, ri4) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(u_par, ri4) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(t_perp, ri4) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(q_perp, ri4) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(q_par, ri4) ! more redundant writes, for backwards compatibility for now: read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(q_par, ri4) read(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) read(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ri4(2,:,:,:,:) = -ri4(2,:,:,:,:) call r2c(q_par, ri4) deallocate (ri4) ! The above replace 3rd and 4th perp moments, not yet used toroidally: ! read(iou,2)((((rps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! read(iou,2)((((rpc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! read(iou,2)((((sps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! read(iou,2)((((spc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) allocate (ri3(2,kdpass, md, nd)) read(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_density, ri3) read(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_p, ri3) read(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_r, ri3) read(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(e_t, ri3) read(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call r2c(phi_ba, ri3) deallocate (ri3) allocate (ri3(2,ld,md,nd)) read(iou,2)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) call r2c(phi_bk, ri3) read(iou,2)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) read(iou,2)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) call r2c(e_denk, ri3) deallocate (ri3) allocate (ri2(2,md,nd)) read(iou,2)((ri2(1,m,n),m=1,md),n=1,nd) read(iou,2)((ri2(2,m,n),m=1,md),n=1,nd) call r2c(phiav_old, ri2) deallocate (ri2) read(iou,2)(timo(i),i=nt1,nt2) !time history plot array read(iou,2)((wpfx(i,n),i=nt1,nt2),n=1,n1) !ion thermal flux read(iou,2)((fluxi(i,n),i=nt1,nt2),n=1,n1) !ion particle flux read(iou,2)(fluxe(i),i=nt1,nt2) !electron particle flux read(iou,2)(qfluxe(i),i=nt1,nt2) !electron thermal flux read(iou,2)(gamx(i),i=nt1,nt2) !total growth rate? read(iou,2)(((mgamx(m,n,i),m=1,md),n=1,nd),i=nt1,nt2) !gamma vs. mode ! read(iou,2)(wrhx(i),i=nt1,nt2) ! read(iou,2)(wrhy(i),i=nt1,nt2) ! read(iou,2)(wdux(i),i=nt1,nt2) ! read(iou,2)(wu2x(i),i=nt1,nt2) ! read(iou,2)(wp2x(i),i=nt1,nt2) ! read(iou,2)(wv2x(i),i=nt1,nt2) ! read(iou,2)(wlux(i),i=nt1,nt2) ! read(iou,2)(wdpx(i),i=nt1,nt2) ! read(iou,2)(wdvx(i),i=nt1,nt2) read(iou,2)(wenx(i),i=nt1,nt2) ! total energy ! read(iou,2)(wdfx(i),i=nt1,nt2) read(iou,2)(pcerr(i),i=nt1,nt2) ! energy conservation read(iou,2)(((wtif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode read(iou,2)(((wkif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !kin.en. vs. mode read(iou,2)(((wpif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !p^2 vs mode read(iou,2)(((psp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !phi*phi vs mode read(iou,2)(((wenrk(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode read(iou,2)(((dketa(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !eta drive v mode read(iou,2)(((dkpar(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !par damp. v mode read(iou,2)(((dktor(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor drive v mode read(iou,2)(((dktdp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor damp. v mode read(iou,2)(((wkups(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !wupsi v mode read(iou,2)(((dktot(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tot. drive v mode read(iou,2)((cx(i,n),i=nt1,nt2),n=1,4*nd) !phi(x) not aved. in y read(iou,2)((cy(i,m),i=nt1,nt2),m=1,4*md) !phi(y) not aved. in x read(iou,2)((cz(i,n),i=nt1,nt2),n=1,ld) !phi(z) not aved. in x,y read(iou,2)((czc(i,n),i=nt1,nt2),n=1,ld) !phi(0)phi(z) aved. in x,y read(iou,2)((czcn(i,n),i=nt1,nt2),n=1,ld) !phi(z)^2 aved. in x,y read(iou,2)((czn(i,n),i=nt1,nt2),n=1,ld) !ne(0)ne(z) aved. in x,y read(iou,2)((cznn(i,n),i=nt1,nt2),n=1,ld) !ne(z)^2 aved. in x,y read(iou,2)(ct(i),i=nt1,nt2) !phi(0)phi(t) aved. in x,y read(iou,1) ninterv read(iou,2)((phi00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. potential read(iou,2)((den00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. density read(iou,2)((upar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. velocity read(iou,2)((tpar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_par read(iou,2)((tperp00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_perp read(iou,2)((drt(id,i),id=1,14),i=nt1,nt2) !drive terms read(iou,2)(phirms(i),i=nt1,nt2) ! rms phi ! read(iou,2)(((wvif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! read(iou,2)(((wplx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) read(iou,2)(wakx(i),i=nt1,nt2) ! mode width in x read(iou,2)(waky(i),i=nt1,nt2) ! mode width in y read(iou,2)(wxsp(i),i=nt1,nt2) ! mode spread read(iou,2)(utor(i),i=nt1,nt2) ! toroidal flow read(iou,2)(upol(i),i=nt1,nt2) ! poloidal flow read(iou,2)(((grtmx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! amplitude allocate (ri3(2, ntf1:ntf2, md, nd)) read(iou,2)(((ri3(1,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) ! Phi(t) read(iou,2)(((ri3(2,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) call r2c(utim(ntf1:ntf2,:,:), ri3(:,ntf1:ntf2,:,:)) deallocate (ri3) read(iou,2)(timf(i),i=ntf1,ntf2) ! time base for freq(t) and phi(t) read(iou,2)avgflag read(iou,3)date read(iou,4)note ! ! Write multi-species stuff at the end for backwards compatibility: ! read(iou,2)(rmass(n),n=1,nspecies) read(iou,2)(charge(n),n=1,nspecies) read(iou,2)(n_I(n),n=1,nspecies) read(iou,2)(Ln(n),n=1,nspecies) read(iou,2)(eta(n),n=1,nspecies) read(iou,2)(eta_par(n),n=1,nspecies) read(iou,2)(tau(n),n=1,nspecies) read(iou,2)(((phisq(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) read(iou,2)((((fluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) read(iou,2)((((qfluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) read(iou,2)(((fluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) read(iou,2)(((qfluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! more new variables tacked on at the end for upward compabitility: read(iou,1) icrit,ntimf read(iou,2)(time(m),m=1,md) read(iou,2)(dt(m),m=1,md) read(iou,2)(((rkperp2(l,m,n),l=1,ld),m=1,md),n=1,nd) read(iou,2) ((rdiff(m,n),m=1,md),n=1,nd) read(iou,2) chi_int read(iou,2) malias, nalias 1 format(1x,3i10) 2 format(1x,3e22.14) 3 format(1x,a30) 4 format(1x,a80) end subroutine wread subroutine wpunch(iou,ntim,ntimr,tim,dt1) ! ! Write/read results to/from the *.res file. ! The format of the *.res file is essentially given by the ! read/writes in this subroutine, and must be compatible between ! output subroutine wpunch.f and the input subroutine wread.f ! ! (c) Copyright 1991 to 1995 by Michael A. Beer, William D. Dorland, ! and Gregory W. Hammett. ALL RIGHTS RESERVED. ! ! subroutine arguments integer iou ! on entry, fortran i/o unit number of the *.res file integer ntim ! on entry, the index of last time point previously used. ! on exit, the index of the last time point now used. integer ntimr ! the same as for ntim, but for a copy of phi(t) at ! one point with potentially higher time resolution ! for frequency calculation. ! output subroutine arguments: real tim ! final time when dt is the same for all ky's. real dt1 ! time step when dt is the same for all ky's. ! local variables: integer i,id,l,m,n,n1,nt1,nt2,ntf1,ntf2 real xxx ! defaults values for backwards compatibility: xxx=0.0 tim=time(md) dt1=dt(md) ! start reading from or writing to the *.res file: write(iou,1)ld,ldb,kd,nmin,nmax,nd,lin,kdpass write(iou,1)nstp,ne,nfreq,ntotal,md,nspecies,ikx,iperiod write(iou,1)nparmom,nperpmom,nemom,iodd,iflr,iphi00,ifilter,iexp,igradon write(iou,1)(mmin(n),n=1,nd) write(iou,1)(mmax(n),n=1,nd) write(iou,2)tim,shr,qsf,epsn,eps,epse,nueeff,alpha, & etai,nuii,rmu1,xxx,xxx,tiovte,vy,vz,rmime,etae write(iou,2)dt0,dt1,x0,y0,z0,xp n1=nspecies nt1=ntim+1 nt2=ntim+ne ntim=ntim+ne ntf1=ntimr+1 ntf2=ntimr+ntotal ntimr=ntimr+ntotal ! eventually, put the following arrays in a more logical order: ! For now, keep backward compability. Put new arrays at end: allocate (ri4(2,ld,md,nd,n1)) call c2r(density, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) allocate (ri3(2,ld,md,nd)) call c2r(potential, ri3) write(iou,2)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) write(iou,2)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) deallocate (ri3) call c2r(t_par, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! redundant write, just to keep other arrays in the same position ! in the file for backwards compatability, for now: call c2r(u_par, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(u_par, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(t_perp, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(q_perp, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(q_par, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! more redundant writes, for backwards compatibility for now: call c2r(q_par, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) call c2r(q_par, ri4) write(iou,2)((((ri4(2,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) write(iou,2)((((ri4(1,l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) deallocate (ri4) ! The above replace 3rd and 4th perp moments, not yet used toroidally: ! write(iou,2)((((rps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! write(iou,2)((((rpc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! write(iou,2)((((sps(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) ! write(iou,2)((((spc(l,m,n,i),l=1,ld),m=1,md),n=1,nd),i=1,n1) allocate (ri3(2,kdpass, md, nd)) call c2r(e_density, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(e_p, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(e_r, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(e_t, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) call c2r(phi_ba, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,kdpass),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,kdpass),m=1,md),n=1,nd) deallocate (ri3) allocate (ri3(2,ld,md,nd)) call c2r(phi_bk, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) call c2r(e_denk, ri3) write(iou,2)(((ri3(1,l,m,n),l=1,ld),m=1,md),n=1,nd) write(iou,2)(((ri3(2,l,m,n),l=1,ld),m=1,md),n=1,nd) deallocate (ri3) allocate (ri2(2,md, nd)) call c2r(phiav_old, ri2) write(iou,2)((ri2(1,m,n),m=1,md),n=1,nd) write(iou,2)((ri2(2,m,n),m=1,md),n=1,nd) deallocate (ri2) write(iou,2)(timo(i),i=nt1,nt2) !time history plot array write(iou,2)((wpfx(i,n),i=nt1,nt2),n=1,n1) !ion thermal flux write(iou,2)((fluxi(i,n),i=nt1,nt2),n=1,n1) !ion particle flux write(iou,2)(fluxe(i),i=nt1,nt2) !electron particle flux write(iou,2)(qfluxe(i),i=nt1,nt2) !electron thermal flux write(iou,2)(gamx(i),i=nt1,nt2) !total growth rate? write(iou,2)(((mgamx(m,n,i),m=1,md),n=1,nd),i=nt1,nt2) !gamma vs. mode ! write(iou,2)(wrhx(i),i=nt1,nt2) ! write(iou,2)(wrhy(i),i=nt1,nt2) ! write(iou,2)(wdux(i),i=nt1,nt2) ! write(iou,2)(wu2x(i),i=nt1,nt2) ! write(iou,2)(wp2x(i),i=nt1,nt2) ! write(iou,2)(wv2x(i),i=nt1,nt2) ! write(iou,2)(wlux(i),i=nt1,nt2) ! write(iou,2)(wdpx(i),i=nt1,nt2) ! write(iou,2)(wdvx(i),i=nt1,nt2) write(iou,2)(wenx(i),i=nt1,nt2) ! total energy ! write(iou,2)(wdfx(i),i=nt1,nt2) write(iou,2)(pcerr(i),i=nt1,nt2) ! energy conservation write(iou,2)(((wtif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode write(iou,2)(((wkif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !kin.en. vs. mode write(iou,2)(((wpif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !p^2 vs mode write(iou,2)(((psp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !phi*phi vs mode write(iou,2)(((wenrk(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !energy vs. mode write(iou,2)(((dketa(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !eta drive v mode write(iou,2)(((dkpar(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !par damp. v mode write(iou,2)(((dktor(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor drive v mode write(iou,2)(((dktdp(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tor damp. v mode write(iou,2)(((wkups(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !wupsi v mode write(iou,2)(((dktot(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) !tot. drive v mode write(iou,2)((cx(i,n),i=nt1,nt2),n=1,4*nd) !phi(x) not aved. in y write(iou,2)((cy(i,m),i=nt1,nt2),m=1,4*md) !phi(y) not aved. in x write(iou,2)((cz(i,n),i=nt1,nt2),n=1,ld) !phi(z) not aved. in x,y write(iou,2)((czc(i,n),i=nt1,nt2),n=1,ld) !phi(0)phi(z) aved. in x,y write(iou,2)((czcn(i,n),i=nt1,nt2),n=1,ld) !phi(z)^2 aved. in x,y write(iou,2)((czn(i,n),i=nt1,nt2),n=1,ld) !ne(0)ne(z) aved. in x,y write(iou,2)((cznn(i,n),i=nt1,nt2),n=1,ld) !ne(z)^2 aved. in x,y write(iou,2)(ct(i),i=nt1,nt2) !phi(0)phi(t) aved. in x,y write(iou,1) ninterv write(iou,2)((phi00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. potential write(iou,2)((den00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. density write(iou,2)((upar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. velocity write(iou,2)((tpar00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_par write(iou,2)((tperp00(i,n),i=nt1,nt2),n=1,nd) !flux surf. ave. t_perp write(iou,2)((drt(id,i),id=1,14),i=nt1,nt2) !drive terms write(iou,2)(phirms(i),i=nt1,nt2) ! rms phi ! write(iou,2)(((wvif(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! write(iou,2)(((wplx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) write(iou,2)(wakx(i),i=nt1,nt2) ! mode width in x write(iou,2)(waky(i),i=nt1,nt2) ! mode width in y write(iou,2)(wxsp(i),i=nt1,nt2) ! mode spread write(iou,2)(utor(i),i=nt1,nt2) ! toroidal flow write(iou,2)(upol(i),i=nt1,nt2) ! poloidal flow write(iou,2)(((grtmx(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! amplitude allocate (ri3(2,ntf1:ntf2,md,nd)) call c2r(utim(ntf1:ntf2,:,:), ri3(:,ntf1:ntf2,:,:)) write(iou,2)(((ri3(1,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) ! Phi(t) write(iou,2)(((ri3(2,i,m,n),i=ntf1,ntf2),m=1,md),n=1,nd) deallocate (ri3) write(iou,2)(timf(i),i=ntf1,ntf2) ! time base for freq(t) and phi(t) write(iou,2)avgflag write(iou,3)date write(iou,4)note ! ! Write multi-species stuff at the end for backwards compatibility: ! write(iou,2)(rmass(n),n=1,nspecies) write(iou,2)(charge(n),n=1,nspecies) write(iou,2)(n_I(n),n=1,nspecies) write(iou,2)(Ln(n),n=1,nspecies) write(iou,2)(eta(n),n=1,nspecies) write(iou,2)(eta_par(n),n=1,nspecies) write(iou,2)(tau(n),n=1,nspecies) write(iou,2)(((phisq(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) write(iou,2)((((fluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) write(iou,2)((((qfluximn(i,m,n,l),i=nt1,nt2),m=1,md),n=1,nd),l=1,nspecies) write(iou,2)(((fluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) write(iou,2)(((qfluxemn(i,m,n),i=nt1,nt2),m=1,md),n=1,nd) ! more new variables tacked on at the end for upward compabitility: write(iou,1) icrit,ntimf write(iou,2)(time(m),m=1,md) write(iou,2)(dt(m),m=1,md) write(iou,2)(((rkperp2(l,m,n),l=1,ld),m=1,md),n=1,nd) write(iou,2) ((rdiff(m,n),m=1,md),n=1,nd) write(iou,2) chi_int write(iou,1) malias, nalias 1 format(1x,3i10) 2 format(1x,3e22.14) 3 format(1x,a30) 4 format(1x,a80) end subroutine wpunch end module io