Subroutine callspecgen(Annode,x,y,Eg,generation)

This subroutine is activated by $usespecgen. Below is an example. Users are free to change the code except the global variables and subroutine interface. This subroutine is saved in libspecgen.f90

module specgenmod
use nrtype
integer(i8b) :: specgen_Nregion
integer(i8b),allocatable :: specgen_type(:)
real(dp),allocatable :: specgen_xc(:),specgen_sgx(:) 
real(dp),allocatable :: specgen_yc(:),specgen_sgy(:)
real(dp),allocatable :: specgen_A(:),specgen_B(:)
logical :: ifspecgen=.false.
end module

subroutine callspecgen(Annode,x,y,Eg,generation)
use nrtype
use global2d
use specgenmod
implicit none
integer(i8b), intent(in) :: Annode
real(dp) , intent(inout) :: x(Annode),y(Annode),eg(Annode),generation(Annode)
real(dp)  :: totalgen
integer :: i,j
  do i=1,Annode
     totalgen=0.0_dp
     do j=1,specgen_Nregion
       if(specgen_type(j)==1) then
         totalgen=totalgen+specgen_A(j)*exp(-(x(i)-specgen_xc(j))**2/(2*specgen_sgx(j)**2)) &
             *specgen_B(j) *exp(-(y(i)-specgen_yc(j))**2/(2*specgen_sgy(j)**2))
        ! print*,totalgen,j,specgen_xc(j),specgen_xy(j),specgen_xgy(j),specgen_sgx(j)
       elseif(specgen_type(j)==2) then
         totalgen=totalgen+specgen_A(j)*exp(-(x(i)-specgen_xc(j))**2/(specgen_sgx(j)**2)) &
             *specgen_B(j)*((y(i)-specgen_yc(j))**2)*exp(-abs(y(i)-specgen_yc(j))/(specgen_sgy(j)))
       elseif(specgen_type(j)==3) then
         totalgen=totalgen+specgen_A(j)*exp(-(y(i)-specgen_yc(j))**2/(specgen_sgy(j)**2)) &
             *specgen_B(j)*((x(i)-specgen_xc(j))**2)*exp(-abs(x(i)-specgen_xc(j))/(specgen_sgx(j)))
       end if
     end do
      generation(i)=generation(i)+totalgen      
  end do
end subroutine