- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi Friends
I have a module which contains two subroutines in which (first subroutine) i use the allocatable attribute to allocate the space for arrays then do some computations and deallocate them at the end. My first subroutine is called more than 1000 times. My question (problem) is when i initilize the allocatable arrays to zero i get a different result and when i don't initilize them to zero i get a different results?.
First of all is it right to initilize the allocatable arrays to zero?. I have seen some threads but iam not able to come to clear conclusion.
I am compiling it in ifort version 12.1.3 Any help is highly appreciated. Below is my code for reference
subroutine rho_z
!===============
implicit none
integer :: i,actual_particle,nrot
integer :: max_pos
real :: ord_xx,ord_yy,ord_zz
real :: ord_xy,ord_xz,ord_yz
integer, dimension(:), allocatable :: index,jndex,counter
real, dimension(:), allocatable :: value_N,value_O,ord_z
real, dimension(:), allocatable :: local_order,help_dip
real, dimension(3,3) :: Q,eigen_vector_glob
real, dimension( 3) :: eigen_value_glob
allocate(index (N_particle)); allocate(local_order (N_particle))
allocate(value_N (N_particle)); allocate(value_O (N_particle))
allocate(ord_z (N_particle))
allocate(jndex (maxchannel_rho));allocate(help_dip(maxchannel_rho))
allocate(counter (maxchannel_rho))
index=0; jndex=0; counter=0 ! Initilizing the allocatable arrays to zero
value_N=0.;value_O=0.;ord_z=0. ! Initilizing the allocatable arrays to zero
local_order=0.;help_dip=0. ! Initilizing the allocatable arrays to zero
Q=0.;eigen_vector_glob=0.
eigen_value_glob=0.
ord_xx=0.0; ord_yy=0.0; ord_zz=0.0
ord_xy=0.0; ord_xz=0.0; ord_yz=0.0
nrot=0
do i=1,N_particle
actual_particle=living_particles%item(i)
index (i)=int((molecule(actual_particle)%z+0.5)/delta_rho)+1
ord_z (i)=molecule(actual_particle)%u_z
value_N(i)=1.0
value_O(i)=0.5*(3.*molecule(actual_particle)%u_z&
*molecule(actual_particle)%u_z -1.)
end do
call av_add_hist(local_density,N_particle,index,value_N)
call av_add_hist(S_local ,N_particle,index,value_O)
call av_add_hist(local_order_z,N_particle,index, ord_z) !Binning the histogram order wrt wall normal
do i=1,N_particle
actual_particle=living_particles%item(i)
ord_xx=ord_xx+0.5*(3.*molecule(actual_particle)%u_x*molecule(actual_particle)%u_x -1.)
ord_yy=ord_yy+0.5*(3.*molecule(actual_particle)%u_y*molecule(actual_particle)%u_y -1.)
ord_zz=ord_zz+0.5*(3.*molecule(actual_particle)%u_z*molecule(actual_particle)%u_z -1.)
ord_xy=ord_xy+1.5* molecule(actual_particle)%u_x*molecule(actual_particle)%u_y
ord_xz=ord_xz+1.5* molecule(actual_particle)%u_x*molecule(actual_particle)%u_z
ord_yz=ord_yz+1.5* molecule(actual_particle)%u_y*molecule(actual_particle)%u_z
end do
ord_xx=ord_xx/N_particle
ord_yy=ord_yy/N_particle
ord_zz=ord_zz/N_particle
ord_xy=ord_xy/N_particle
ord_xz=ord_xz/N_particle
ord_yz=ord_yz/N_particle
Q(1,1)=ord_xx
Q(2,2)=ord_yy
Q(3,3)=ord_zz
Q(1,2)=ord_xy
Q(2,1)=ord_xy
Q(1,3)=ord_xz
Q(3,1)=ord_xz
Q(2,3)=ord_yz
Q(3,2)=ord_yz
call Jacobi(Q,3,eigen_value_glob,eigen_vector_glob,nrot)
max_pos = maxloc(eigen_value_glob,1)
do i=1,N_particle
actual_particle=living_particles%item(i)
index (i)=int((molecule(actual_particle)%z+0.5)/delta_rho)+1
local_order (i)=molecule(actual_particle)%u_x*eigen_vector_glob(1,max_pos)+&
molecule(actual_particle)%u_y*eigen_vector_glob(2,max_pos)+&
molecule(actual_particle)%u_z*eigen_vector_glob(3,max_pos)
help_dip (index(i))=help_dip(index(i))+local_order(i)
counter (index(i))=counter(index(i)) + 1 !Counter for no.of particles in a bin
!write(*,*)help_dip(i),counter(i)
end do
do i=1,maxchannel_rho
jndex(i)=i
if(counter(i).gt.0) then
help_dip(i)=abs(help_dip(i))/float(counter(i))
else
help_dip(i)=0.
endif
!write(6,*)counter(i)
end do
call av_add_hist(local_dipol_or,maxchannel_rho,jndex,help_dip) !Binning the histogram local dipolar parameter
deallocate(index) ;deallocate(jndex)
deallocate(value_N) ;deallocate(value_O)
deallocate(help_dip);deallocate(ord_z)
deallocate(counter);deallocate(local_order)
return
end subroutine rho_z
Link Copied
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
First of all is it right to initilize the allocatable arrays to zero?. I have seen some threads but iam not able to come to clear conclusion.
Arrays as well as scalar variables, whether allocated on the heap or the stack, have to be initialized before they are used. The programmer is responsible for doing the initialization, which can be done by (i) executable statements that assign values to the variables (assignment statements, READ statements, ...), (ii) initialization by DATA statements or (iii) default initialization of user defined types.
Zero is often not a proper initialization for a newly allocated variable. For example, character variables and logical variables cannot be initialized to zero. Likewise as to the speed of light, the age of the Earth, the number of atoms in the universe, etc. The proper initialization value is quite domain dependent, and you need to judge what is proper for your application.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
@mecej4 Thank you very much for your clarification. To me it, still looks like zero is the proper initilization for my application because i am going to keep on adding numbers between 0 to 1 in the repective array elements. Kindly, let me know if you have any comments.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Just a comment on mecej4's reply - DATA statements can't be used with allocatable arrays.
vinoth, only you can determine whether 0 is the correct initial value. It seems that in this case it is. But the fundamental issue is that you are responsible for making sure that every variable is properly initialized before it is referenced.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Thank you Steve for your inputs.
- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page