Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Ptc/10 #522

Merged
merged 4 commits into from
Sep 28, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
201 changes: 144 additions & 57 deletions forest/code/Ci_tpsa.f90

Large diffs are not rendered by default.

16 changes: 10 additions & 6 deletions forest/code/Se_status.f90
Original file line number Diff line number Diff line change
Expand Up @@ -144,7 +144,7 @@ module S_status
private track_TREE_G_complexr,track_TREE_G_complexp,track_TREE_probe_complexr,track_TREE_probe_complexp_new
real(dp),TARGET ::INITIAL_CHARGE=1
logical :: mcmillan=.false.
real(dp) :: radfac=1 ! to fudge radiation (lower it)
real(dp) :: radfac=1,flucfac=1 ! to fudge radiation (lower it)
TYPE B_CYL
integer firsttime
integer, POINTER :: nmul,n_mono !,nmul_e,n_mono_e
Expand Down Expand Up @@ -264,7 +264,7 @@ real(dp) function cflucf(p)
! if(junk_e) then
! cflucf=cfluc*p%p0c**5
! else
cflucf=cfluc0*twopii/p%GAMMA0I**5/p%MASS**2
cflucf=flucfac*cfluc0*twopii/p%GAMMA0I**5/p%MASS**2
! endif
end function cflucf

Expand Down Expand Up @@ -1529,6 +1529,8 @@ subroutine S_init(STATE,NO1,NP1,pack,ND2,NPARA,number_of_clocks)
INTEGER,optional :: ND2,NPARA,number_of_clocks
INTEGER ND2l,NPARAl,n_acc,no1c,nv,i
LOGICAL(lp) package


do_damping=.false.
do_spin=.false.
if(state%radiation) do_damping=.true.
Expand Down Expand Up @@ -1589,7 +1591,7 @@ subroutine S_init(STATE,NO1,NP1,pack,ND2,NPARA,number_of_clocks)
! endif
! write(6,*) NO1,ND1,NP1,NDEL,NDPT1
!pause 678

CALL INIT(NO1,ND1,NP1+NDEL+2*n_acc,NDPT1,PACKAGE)
nv=2*nd1+NP1+NDEL+2*n_acc

Expand All @@ -1611,12 +1613,14 @@ subroutine S_init(STATE,NO1,NP1,pack,ND2,NPARA,number_of_clocks)
ND1=ND1+n_acc
if(use_complex_in_ptc) call c_init(NO1c,nd1,np1+ndel,ndpt1,n_acc,ptc=my_false) ! PTC false because we will not use the real FPP for acc modulation
n_rf=n_acc

previous_newtpsa=newtpsa
END subroutine S_init

subroutine kill_map_cp()
implicit none

logical present_newtpsa
present_newtpsa=newtpsa
newtpsa=previous_newtpsa
if(associated(dz_8)) then
call kill(dz_8)
deallocate(dz_8)
Expand All @@ -1635,7 +1639,7 @@ subroutine kill_map_cp()
deallocate(dz_c)
nullify(dz_c)
endif

newtpsa=present_newtpsa
end subroutine kill_map_cp


Expand Down
68 changes: 34 additions & 34 deletions forest/code/Si_def_element.f90
Original file line number Diff line number Diff line change
Expand Up @@ -2031,68 +2031,68 @@ SUBROUTINE transfer_ANBN(EL,ELP,VR,DVR,VP,DVP,T)

if(EL%KIND==kind1) return

if(associated(EL%ramp)) then

if(EL%KIND/=kind15) then
do n=1,EL%P%NMUL
EL%BN(N)= EL%ramp%table(0)%bn(n)
EL%AN(N)= EL%ramp%table(0)%an(n)
ELP%BN(N)= ELP%ramp%table(0)%bn(n)
ELP%AN(N)= ELP%ramp%table(0)%an(n)
enddo
else
EL%VOLT=EL%ramp%table(0)%bn(1)*COS(twopi*EL%ramp%table(0)%an(1)*T/clight+EL%ramp%table(0)%bn(2))+EL%ramp%table(0)%an(2)
ELP%VOLT=EL%ramp%table(0)%bn(1)*COS(twopi*EL%ramp%table(0)%an(1)*T/clight+EL%ramp%table(0)%bn(2))+EL%ramp%table(0)%an(2)
write(6,*) " volt ",el%volt,EL%ramp%table(0)%bn(1)
endif

if(EL%ramp%table(0)%b_t/=0.0_dp) then
if(EL%parent_fibre%PATCH%TIME==0) EL%parent_fibre%PATCH%TIME=2
if(EL%parent_fibre%PATCH%TIME==1) EL%parent_fibre%PATCH%TIME=3
EL%parent_fibre%PATCH%b_T=EL%ramp%table(0)%b_t
else
if(EL%parent_fibre%PATCH%TIME==2) EL%parent_fibre%PATCH%TIME=0
if(EL%parent_fibre%PATCH%TIME==3) EL%parent_fibre%PATCH%TIME=1
EL%parent_fibre%PATCH%b_T=0.0_dp
endif
! if(associated(EL%ramp)) then
!
! if(EL%KIND/=kind15) then
! do n=1,EL%P%NMUL
! EL%BN(N)= EL%ramp%table(0)%bn(n)
! EL%AN(N)= EL%ramp%table(0)%an(n)
! ELP%BN(N)= ELP%ramp%table(0)%bn(n)
! ELP%AN(N)= ELP%ramp%table(0)%an(n)
! enddo
! else
! EL%VOLT=EL%ramp%table(0)%bn(1)*COS(twopi*EL%ramp%table(0)%an(1)*T/clight+EL%ramp%table(0)%bn(2))+EL%ramp%table(0)%an(2)
! ELP%VOLT=EL%ramp%table(0)%bn(1)*COS(twopi*EL%ramp%table(0)%an(1)*T/clight+EL%ramp%table(0)%bn(2))+EL%ramp%table(0)%an(2)
! write(6,*) " volt ",el%volt,EL%ramp%table(0)%bn(1)
! endif
!
! if(EL%ramp%table(0)%b_t/=0.0_dp) then
! if(EL%parent_fibre%PATCH%TIME==0) EL%parent_fibre%PATCH%TIME=2
! if(EL%parent_fibre%PATCH%TIME==1) EL%parent_fibre%PATCH%TIME=3
! EL%parent_fibre%PATCH%b_T=EL%ramp%table(0)%b_t
! else
! if(EL%parent_fibre%PATCH%TIME==2) EL%parent_fibre%PATCH%TIME=0
! if(EL%parent_fibre%PATCH%TIME==3) EL%parent_fibre%PATCH%TIME=1
! EL%parent_fibre%PATCH%b_T=0.0_dp
! endif

else
! else

IF(EL%P%NMUL>=1.and.associated(EL%D_BN)) THEN
if(present(VR))then
do n=1,EL%P%NMUL
EL%BN(N)= vR*EL%D0_BN(N)+DVR*EL%D_BN(N)
EL%AN(N)= vR*EL%D0_AN(N)+DVR*EL%D_AN(N)
ELP%BN(N)= vR*EL%D0_BN(N)+DVR*EL%D_BN(N)
ELP%AN(N)= vR*EL%D0_AN(N)+DVR*EL%D_AN(N)
ELP%BN(N)= vR*ELp%D0_BN(N)+DVR*ELp%D_BN(N)
ELP%AN(N)= vR*ELp%D0_AN(N)+DVR*ELp%D_AN(N)
enddo
else
do n=1,EL%P%NMUL
EL%BN(N)= vp*EL%D0_BN(N)+DVp*EL%D_BN(N)
EL%AN(N)= vp*EL%D0_AN(N)+DVp*EL%D_AN(N)
ELP%BN(N)= vp*EL%D0_BN(N)+DVp*EL%D_BN(N)
ELP%AN(N)= vp*EL%D0_AN(N)+DVp*EL%D_AN(N)
ELP%BN(N)= vp*ELp%D0_BN(N)+DVp*ELp%D_BN(N)
ELP%AN(N)= vp*ELp%D0_AN(N)+DVp*ELp%D_AN(N)
enddo
endif


endif
! endif
if(associated(el%volt)) then
if(present(VR))then
EL%volt= vR*EL%D0_Volt+DVR*EL%D_Volt
ELP%volt= vR*EL%D0_Volt+DVR*EL%D_Volt
ELP%volt= vR*ELp%D0_Volt+DVR*ELp%D_Volt
else
EL%volt= vp*EL%D0_Volt+DVp*EL%D_Volt
ELP%volt= vp*EL%D0_Volt+DVp*EL%D_Volt
ELP%volt= vp*ELp%D0_Volt+DVp*ELp%D_Volt
endif
endif
if(associated(el%phas)) then
if(present(VR))then
EL%phas= vR*EL%D0_phas+DVR*EL%D_phas
ELP%phas= vR*EL%D0_phas+DVR*EL%D_phas
ELP%phas= vR*ELp%D0_phas+DVR*ELp%D_phas
else
EL%phas= vp*EL%D0_phas+DVp*EL%D_phas
ELP%phas= vp*EL%D0_phas+DVp*EL%D_phas
ELP%phas= vp*ELp%D0_phas+DVp*ELp%D_phas
endif
endif

Expand Down
4 changes: 4 additions & 0 deletions forest/code/Sr_spin.f90
Original file line number Diff line number Diff line change
Expand Up @@ -1572,8 +1572,10 @@ SUBROUTINE TRACK_NODE_FLAG_probe_QUAP(C,XS,K)
! endif

IF(K%MODULATION.and.xs%nac/=0) then
knob=.true.
if(c%parent_fibre%mag%slow_ac/=0) CALL MODULATE(C,XS,K) !modulate
CALL TRACK_MODULATION(C,XS,K) !modulate
knob=.false.
ENDIF !modulate


Expand Down Expand Up @@ -1939,8 +1941,10 @@ SUBROUTINE TRACK_NODE_FLAG_probe_P(C,XS,K)
! endif

IF(K%MODULATION.and.xs%nac/=0) then
knob=.true.
if(c%parent_fibre%mag%slow_ac/=0) CALL MODULATE(C,XS,K) !modulate
CALL TRACK_MODULATION(C,XS,K) !modulate
knob=.false.
ENDIF !modulate


Expand Down
12 changes: 7 additions & 5 deletions forest/code/Su_duan_zhe_map.f90
Original file line number Diff line number Diff line change
Expand Up @@ -1073,18 +1073,19 @@ end SUBROUTINE orthonormalise

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! new zhe tracking !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

SUBROUTINE track_TREE_probe_complex_zhe(T,xs,spin,rad,stoch,linear,slim)
SUBROUTINE track_TREE_probe_complex_zhe(T,xs,spin,rad,stoch,linear,slim,flucfac)
! use da_arrays
IMPLICIT NONE
TYPE(TREE_ELEMENT),target, INTENT(INout) :: T(3)

type(probe) xs
real(dp) x(size_tree),x0(size_tree),s0(3,3),r(3,3),dx6,beta,q(3),p(3),qg(3),qf(3)
real(dp) normb,norm,x0_begin(size_tree),xr(6),normbb
real(dp) normb,norm,x0_begin(size_tree),xr(6),normbb,flucfac1
integer i,j,k,ier,is
logical, optional :: spin,stoch,rad,linear,slim
logical spin0,stoch0,rad0,doit,as_is0,slim0
integer no1
real(dp), optional :: flucfac
type(quaternion) qu
as_is0=t(1)%usenonsymp
spin0=.true.
Expand All @@ -1099,7 +1100,8 @@ SUBROUTINE track_TREE_probe_complex_zhe(T,xs,spin,rad,stoch,linear,slim)
if(present(spin)) spin0=spin
if(present(stoch)) stoch0=stoch
if(present(rad)) rad0=rad

flucfac1=1
if(present(flucfac)) flucfac1=flucfac
if(as_is0) rad0=.true.
doit=rad0.or.stoch0
x=0
Expand All @@ -1117,11 +1119,11 @@ SUBROUTINE track_TREE_probe_complex_zhe(T,xs,spin,rad,stoch,linear,slim)
xr=0.0_dp
if(use_gaussian_zhe) then
do i=1,6
xr(i)=GRNF_zhe_gaussian()*t(2)%fix0(i)
xr(i)=flucfac1*GRNF_zhe_gaussian()*t(2)%fix0(i)
enddo
else
do i=1,6
xr(i)=GRNF_zhe()*t(2)%fix0(i)
xr(i)=flucfac1*GRNF_zhe()*t(2)%fix0(i)
enddo
endif
xr =matmul(t(2)%rad,xr)
Expand Down
13 changes: 12 additions & 1 deletion forest/code/a_scratch_size.f90
Original file line number Diff line number Diff line change
Expand Up @@ -226,6 +226,8 @@ module precision_constants
character(18) :: format8="(8(1x,g23.16,1x))"
character(18) :: format9="(9(1x,g23.16,1x))"
character(19) :: format10="(10(1x,g23.16,1x))"
character(19) :: format11="(11(1x,g23.16,1x))"
character(19) :: format12="(12(1x,g23.16,1x))"
! logical(lp) :: fixed_found
! lielib_print(1)=1 lieinit prints info
! lielib_print(2)=1 expflo warning if no convergence
Expand All @@ -246,7 +248,16 @@ module precision_constants
! lielib_print(17)=1 print magnets with excessive cutting
integer , target :: old_integrator =1 ! before making spin high order
character*255 :: preffile="pref.txt"

logical :: newtpsa=.false., assume_da_map = .true.,previous_newtpsa=.false.,use_np=.true.,newspin=.true.
integer :: with_para=2
integer :: nphere=0
integer(4), pointer :: inds(:,:),ind1(:),ind2(:) => null()
integer(4), pointer :: nind1(:),nind2(:) => null()
integer(4) combien, poscombien,ninds
real(dp), pointer:: reel(:),finds(:,:,:) => null(),finds1(:) => null()
complex(dp), pointer:: reelc(:) => null()
integer :: ipric=0,ipric2=0

INTERFACE read
MODULE PROCEDURE read_d
MODULE PROCEDURE read_int
Expand Down
Loading