module nonlinear ! ! (c) Copyright 1991 to 1998 by Michael A. Beer, William D. Dorland, ! P. B. Snyder, Q. P. Liu, and Gregory W. Hammett. ALL RIGHTS RESERVED. ! ! ! This module contains the nonlinear terms. ! implicit none ! ! Nonlinear terms. Could be condensed down to one per equation for ! a small memory cost. ! complex, dimension(:,:,:,:), allocatable :: nl_density1, nl_density2, & nl_upar1, nl_upar2, nl_tpar, nl_qpar, nl_tperp1, nl_tperp2, & nl_tperp3, nl_qperp1, nl_qperp2, nl_qperp3 ! ! Nonlinear finite beta terms ! complex, dimension(:,:,:,:), allocatable :: emnl_density1, & emnl_upar1, emnl_upar2, emnl_tpar, emnl_tperp, & uparnl, tparnl complex, dimension(:,:,:), allocatable :: emnl_n_e1, emnl_n_e2, & emnl_apar, emnl_t_e ! ! Electron nonlinear terms in pitch angle space ! complex, dimension(:,:,:), allocatable :: nl_edensity, nl_ep, nl_er, nl_et contains subroutine nonlin(density1, u_par1, t_par1, q_par1, t_perp1, q_perp1, phi) use itg_data, only: nspecies, nparmom, nperpmom, iflr use nlps, only: nlpsc use gryffin_grid, only: ldb use linear, only: phi_flr1, phi_flr2, phi_u use gryffin_layouts use mp, only: iproc complex, dimension(:,m_low:,n_low:) :: phi complex, dimension(:,m_low:,n_low:,:) :: & density1, u_par1, t_par1, q_par1, t_perp1, q_perp1 integer i ! ! We must calculate {a,b} for a,b=: ! uv,w uv,tpar uv,v uv,T ! uv,q uv,qpar uflr,w uflr,T ! uflr,v uflr2,T uflr2,q ! ! In each case we want to pass the two variables' real and ! imaginary parts to a subroutine which will return the real ! and imaginary parts of {a,b} in (x,ky,kz) space. ! ! ! Usual NL terms unless iflr=6, then do like Ron ! if (iflr /= 6) then do i=1,nspecies call nlpsc(0,phi_flr2(:,:,:,i),t_perp1(:,:,:,i),nl_tperp3(:,:,:,i),ldb) call nlpsc(2,phi_flr1(:,:,:,i),t_perp1(:,:,:,i),nl_density2(:,:,:,i),ldb) call nlpsc(1,phi_flr1(:,:,:,i),density1(:,:,:,i),nl_tperp2(:,:,:,i),ldb) call nlpsc(2,phi_u(:,:,:,i),density1(:,:,:,i),nl_density1(:,:,:,i),ldb) call nlpsc(1,phi_u(:,:,:,i),t_par1(:,:,:,i),nl_tpar(:,:,:,i),ldb) call nlpsc(1,phi_u(:,:,:,i),u_par1(:,:,:,i),nl_upar1(:,:,:,i),ldb) call nlpsc(1,phi_u(:,:,:,i),t_perp1(:,:,:,i),nl_tperp1(:,:,:,i),ldb) if(nparmom == 4) call nlpsc(1,phi_u(:,:,:,i),q_par1(:,:,:,i),nl_qpar(:,:,:,i),ldb) if(nperpmom == 2) then call nlpsc(1,phi_u(:,:,:,i),q_perp1(:,:,:,i),nl_qperp1(:,:,:,i),ldb) call nlpsc(2,phi_flr2(:,:,:,i),q_perp1(:,:,:,i),nl_qperp3(:,:,:,i),ldb) call nlpsc(2,phi_flr1(:,:,:,i),q_perp1(:,:,:,i),nl_upar2(:,:,:,i),ldb) call nlpsc(1,phi_flr1(:,:,:,i),u_par1(:,:,:,i),nl_qperp2(:,:,:,i),ldb) endif enddo else ! ! Like Ron (3+1 assumed) ! do i=1,nspecies call nlpsc(0,phi,density1(:,:,:,i),nl_density1(:,:,:,i),ldb) call nlpsc(1,phi,u_par1(:,:,:,i),nl_upar1(:,:,:,i),ldb) call nlpsc(1,phi,t_par1(:,:,:,i),nl_tpar(:,:,:,i),ldb) call nlpsc(1,phi,t_perp1(:,:,:,i),nl_tperp1(:,:,:,i),ldb) enddo endif end subroutine nonlin subroutine nonline(e_density1, e_p1, e_r1, phi_ba1) use gryffin_grid, only: kd use nlps, only: nlpsc use gryffin_layouts complex, dimension(:,m_low:,n_low:) :: & phi_ba1, e_density1, e_p1, e_r1 call nlpsc(0,phi_ba1,e_density1,nl_edensity,kd) call nlpsc(1,phi_ba1,e_p1,nl_ep,kd) call nlpsc(1,phi_ba1,e_r1,nl_er,kd) end subroutine nonline subroutine nonlin_em(density1, u_par1, t_par1, q_par1, t_perp1, & q_perp1, phi, n_e1, u_e1, tpar_e1, tperp_e1, apar1) use itg_data, only: iflr, nspecies, tiovte use gryffin_grid, only: ldb use fields, only: n_e, u_e, phi_na use gryffin_layouts use nlps, only: nlpsc use linear, only: phi_flr1, phi_flr2, phi_u, J0apar complex, dimension(:,m_low:,n_low:) :: phi complex, dimension(:,m_low:,n_low:,:) :: & density1, u_par1, t_par1, q_par1, t_perp1, q_perp1 complex, dimension(:,m_low:,n_low:) :: & n_e1, u_e1, tpar_e1, tperp_e1, apar1 integer i, l, m, n ! electromagnetic passing electron terms ! (note that terms will subtracted, not added) do n = n_low, n_high do m = m_low, m_high phi_na(:,m,n) = phi(:,m,n) - n_e1(:,m,n)/tiovte do i = 1, nspecies uparnl(:,m,n,i) = density1(:,m,n,i)+t_par1(:,m,n,i) & +phi_u(:,m,n,i) tparnl(:,m,n,i)=2.*u_par1(:,m,n,i)+q_par1(:,m,n,i) enddo enddo enddo call nlpsc(0, phi, n_e1, emnl_n_e1, ldb) call nlpsc(0, apar1, u_e1, emnl_n_e2, ldb) call nlpsc(1, apar1, phi_na, emnl_apar, ldb) call nlpsc(1, apar1, tpar_e1, emnl_t_e, ldb) if (iflr /= 6) then do i=1,nspecies call nlpsc(0,J0apar(:,:,:,i),u_par1(:,:,:,i),emnl_density1(:,:,:,i),ldb) call nlpsc(1,J0apar(:,:,:,i),uparnl(:,:,:,i),emnl_upar1(:,:,:,i),ldb) call nlpsc(1,J0apar(:,:,:,i),tparnl(:,:,:,i),emnl_tpar(:,:,:,i),ldb) call nlpsc(1,J0apar(:,:,:,i),q_perp1(:,:,:,i),emnl_tperp(:,:,:,i),ldb) enddo else ! ! Like Ron (3+1 assumed) ! do i=1,nspecies call nlpsc(0,apar1,u_par1(:,:,:,i),emnl_density1(:,:,:,i),ldb) enddo endif end subroutine nonlin_em subroutine alloc_nonlinear use itg_data, only: epse, beta_e, nspecies use gryffin_layouts use gryffin_grid, only: ld, kd allocate (nl_density1(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_density2(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_upar1(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_upar2(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_tpar(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_qpar(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_tperp1(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_tperp2(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_tperp3(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_qperp1(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_qperp2(ld, m_low:m_alloc, n_low:n_alloc, nspecies), & nl_qperp3(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) nl_density1 = 0.; nl_density2 = 0.; nl_upar1 = 0.; nl_upar2 = 0. nl_tpar = 0.; nl_qpar = 0.; nl_tperp1 = 0.; nl_tperp2 = 0. nl_tperp3 = 0.; nl_qperp1 = 0.; nl_qperp2 = 0.; nl_qperp3 = 0. allocate (emnl_density1(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) allocate (emnl_upar1(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) allocate (emnl_upar2(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) allocate (emnl_tpar(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) allocate (emnl_tperp(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) emnl_density1 = 0.; emnl_upar1 = 0. emnl_upar2 = 0.; emnl_tpar = 0.; emnl_tperp = 0. if(epse > 0.) then allocate (nl_edensity(kd, m_low:m_alloc, n_low:n_alloc), & nl_ep(kd, m_low:m_alloc, n_low:n_alloc), & nl_er(kd, m_low:m_alloc, n_low:n_alloc), & nl_et(kd, m_low:m_alloc, n_low:n_alloc)) nl_edensity = 0.; nl_ep = 0.; nl_er = 0.; nl_et = 0. endif if (beta_e > 0.) then allocate (uparnl(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) allocate (tparnl(ld, m_low:m_alloc, n_low:n_alloc, nspecies)) uparnl = 0.; tparnl = 0. allocate (emnl_n_e1(ld, m_low:m_alloc, n_low:n_alloc)) allocate (emnl_n_e2(ld, m_low:m_alloc, n_low:n_alloc)) allocate (emnl_apar(ld, m_low:m_alloc, n_low:n_alloc)) allocate (emnl_t_e(ld, m_low:m_alloc, n_low:n_alloc)) emnl_n_e1 = 0.; emnl_n_e2 = 0.; emnl_apar = 0.; emnl_t_e = 0. endif end subroutine alloc_nonlinear end module nonlinear