コミット: c4a68e721236d95bb9a398678120b101126669b3 - myprojects-hg-reborn (hg)

コミットメタ情報

• リビジョン: c4a68e721236d95bb9a398678120b101126669b3 (tree)
• 日時: 2009-01-15 23:59:04
• 作者: iselllo
• コミッター: iselllo

ログメッセージ

I wrote a new fortran routine which again calculates molecule absorption for a given system snapshot.

変更サマリ

• add: f2py/absorbe_snapshot.f90

差分

diff -r 56b1b371e4f5 -r c4a68e721236 f2py/absorbe_snapshot.f90
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/f2py/absorbe_snapshot.f90	Thu Jan 15 14:59:04 2009 +0000
@@ -0,0 +1,146 @@
+subroutine absorption_snap_calc(trajectories_arr,N_trajectories, x_arr,y_arr, &
+ z_arr, r1_arr, N_spheres,k, results, old_results)
+
+implicit none
+
+integer ::  N_trajectories, i, j, k, N_spheres ! , count_abso, count_abso_0 ! , step_absorption, flag 
+
+real(kind=8) :: trajectories_arr((3*N_trajectories)), x_arr(N_spheres),y_arr(N_spheres)
+!real(kind=8) :: t_abso_bis
+real(kind=8) :: my_dist, z_arr(N_spheres), surviving_molecules(N_trajectories),r1_arr(N_spheres)
+integer, intent(out) :: results(2,N_trajectories)
+integer  :: old_results(2,N_trajectories)
+
+
+!compile it using:
+! f2py -c --fcompiler=gfortran  absorbe_snapshot.f90 -m --f90flags='-O2' absorbe_snap
+
+
+
+!  count_abso=0
+!  count_abso_0=0
+
+! I will initialize this array within the main code now
+
+!initialize the array with the results
+
+
+
+! do i=1, 2
+
+! do j=1, N_trajectories
+
+! results(i,j)=-1  !I initialize this array with an "impossible" value
+
+! enddo
+! enddo
+
+
+
+
+results=old_results
+
+
+
+!Now check if each single ray touches a circle (first neighbor of our circle)
+
+
+
+trajectory_loop : do i=1, N_trajectories  !go on all the molecule trajectories
+
+!I need to be careful here: I should check only trajectories of molecule which have not been absorbed yet.
+
+
+if (old_results(1,i) .eq. -1) then
+
+!if (k .eq. 20) then
+
+!print * ,"old_results(1,i) == -1 and k", old_results(1,i) == -1, k
+
+!endif
+
+!Now I introduce a condition: I am looking for collisions only for
+!molecules which have NOT already been absorbed. The condition of non-absorption
+!means that surviving_molecules(i) is positive
+
+!if (surviving_molecules(i)>dble(0.0d0)) then
+
+!print *, "Dealing with trajectory number,", i
+
+
+!time_loop : do k=1, NN  !go down in time
+
+!!!! It is important not to swap the loops on time and monomers, otherwise I may end up
+!! seeing a molecule collide with the wrong monomer (in reality it would have been
+!! absorbed at a previous time)
+
+
+
+monomer_loop : do j=1, N_spheres  !check for collisions with a monomer
+
+
+
+
+my_dist=dsqrt(((trajectories_arr((3*i-2))-x_arr(j)))*(trajectories_arr((3*i-2))-x_arr(j))+ &
+  (trajectories_arr((3*i-1))-y_arr(j))*(trajectories_arr((3*i-1))-y_arr(j)) + &
+(trajectories_arr((3*i))-z_arr(j))*(trajectories_arr((3*i))-z_arr(j)))
+
+! if (k .eq. 0) then
+
+! print *, "k, r1_arr(j),x_arr(j),y_arr(j),z_arr(j), and my_dist are", k, r1_arr(j),x_arr(j),y_arr(j),z_arr(j),  my_dist
+
+! print *, "traj((3*i-2)),traj((3*i-1)), traj((3*i)) are, ", trajectories_arr((3*i-2)), &
+! trajectories_arr((3*i-1)), trajectories_arr((3*i))
+
+! endif
+
+if (my_dist .le. r1_arr(j)) then
+
+!count_abso=count_abso+1
+
+ ! if (k .eq. 1) then
+
+!  count_abso_0=count_abso_0+1
+
+!  endif
+
+!print * , "my_dist is, traj,  time, count_abso, count_abso_0 are ", my_dist, i, k, count_abso,count_abso_0 
+
+
+results(1,i)=k  !absorption time (which is now given as an input; I do not have time any longer in this routine)
+results(2,i)=j !id of the absorbing monomer
+
+end if
+
+
+!if (my_dist .le. r1_arr(j)) exit monomer_loop  !leave loop on monomers; the molecule has been absorbed
+!if (my_dist .le. -25.3) exit
+
+
+enddo monomer_loop
+
+
+
+
+
+!if (my_dist .le. r1_arr(j)) exit time_loop  ! I should also leave the loop on time and go back to the trajectories
+
+!if (my_dist .le. -25.3) exit
+
+
+!enddo time_loop
+
+endif
+
+
+enddo trajectory_loop
+
+
+
+
+
+
+return
+
+end subroutine absorption_snap_calc
+