diff options
author | diener <diener@2a26948c-0e4f-0410-aee8-f1d3e353619c> | 2003-09-01 12:39:42 +0000 |
---|---|---|
committer | diener <diener@2a26948c-0e4f-0410-aee8-f1d3e353619c> | 2003-09-01 12:39:42 +0000 |
commit | e54252a77c9efde203c425f84bf4b3bd61d11278 (patch) | |
tree | 98063d923bb0ca9d180d06f24280320ad6f82dbb | |
parent | 77500494f954c4690b6e98f9121b2faab5c48261 (diff) |
Added support for using all the different kinds of initial guesses for the
level set function for all the different level set functions.
git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/EHFinder/trunk@135 2a26948c-0e4f-0410-aee8-f1d3e353619c
-rw-r--r-- | src/EHFinder_Init.F90 | 56 |
1 files changed, 28 insertions, 28 deletions
diff --git a/src/EHFinder_Init.F90 b/src/EHFinder_Init.F90 index 8fc9587..c291376 100644 --- a/src/EHFinder_Init.F90 +++ b/src/EHFinder_Init.F90 @@ -42,10 +42,10 @@ subroutine EHFinder_Init_F(CCTK_ARGUMENTS) deallocate ( reparam_this_level_set ) end if allocate ( reparam_this_level_set(eh_number_level_sets) ) - if ( allocated(reparam_undone) ) then - deallocate ( reparam_undone ) + if ( allocated(re_initialize_undone) ) then + deallocate ( re_initialize_undone ) end if - allocate ( reparam_undone(eh_number_level_sets) ) + allocate ( re_initialize_undone(eh_number_level_sets) ) if ( evolve_generators .gt. 0 ) then @@ -94,16 +94,16 @@ subroutine EHFinder_Init_F(CCTK_ARGUMENTS) ! Set up a sphere of radius initial_rad and translated ! (translate_x,translate_y,translate_z) away from the origin. - f(:,:,:,l) = sqrt( ( x - translate_x )**2 + & - ( y - translate_y )**2 + & - ( z - translate_z )**2 ) - initial_rad(l) + f(:,:,:,l) = sqrt( ( x - translate_x(l) )**2 + & + ( y - translate_y(l) )**2 + & + ( z - translate_z(l) )**2 ) - initial_rad(l) if ( evolve_generators .gt. 0 ) then if ( CCTK_EQUALS( generator_distribution, 'line' ) ) then do i = 1, lsh(1) theta = thetamin + dtheta * ( i + lbnd(1) - 1 ) - xg(i,l) = initial_rad(l) * sin(theta) + translate_x - yg(i,l) = translate_y - zg(i,l) = initial_rad(l) * cos(theta) + translate_z + xg(i,l) = initial_rad(l) * sin(theta) + translate_x(l) + yg(i,l) = translate_y(l) + zg(i,l) = initial_rad(l) * cos(theta) + translate_z(l) end do end if end if @@ -113,12 +113,12 @@ subroutine EHFinder_Init_F(CCTK_ARGUMENTS) if ( CCTK_EQUALS( initial_f(l), 'ellipsoid' ) ) then ! Calculate sines and cosines of the rotation parameters. - cosa = cos(rotation_alpha) - sina = sin(rotation_alpha) - cosb = cos(rotation_beta) - sinb = sin(rotation_beta) - cosc = cos(rotation_gamma) - sinc = sin(rotation_gamma) + cosa = cos(rotation_alpha(l)) + sina = sin(rotation_alpha(l)) + cosb = cos(rotation_beta(l)) + sinb = sin(rotation_beta(l)) + cosc = cos(rotation_gamma(l)) + sinc = sin(rotation_gamma(l)) ! Set up the rotation matrix. The order is alpha around the z-axis, ! beta around the y-axis and finally gamma around the x-axis. @@ -139,13 +139,13 @@ subroutine EHFinder_Init_F(CCTK_ARGUMENTS) do k = 1, nz do j = 1, ny do i = 1, nx - xp(1) = x(i,j,k) - translate_x - xp(2) = y(i,j,k) - translate_y - xp(3) = z(i,j,k) - translate_z + xp(1) = x(i,j,k) - translate_x(l) + xp(2) = y(i,j,k) - translate_y(l) + xp(3) = z(i,j,k) - translate_z(l) xpt = matmul ( txyz, xp ) - f(i,j,k,l) = sqrt( xpt(1)**2 / initial_a**2 + & - xpt(2)**2 / initial_b**2 + & - xpt(3)**2 / initial_c**2) - 1.0 + f(i,j,k,l) = sqrt( xpt(1)**2 / initial_a(l)**2 + & + xpt(2)**2 / initial_b(l)**2 + & + xpt(3)**2 / initial_c(l)**2) - 1.0 end do end do end do @@ -154,11 +154,11 @@ subroutine EHFinder_Init_F(CCTK_ARGUMENTS) if ( CCTK_EQUALS( generator_distribution, 'line' ) ) then do i = 1, lsh(1) theta = thetamin + dtheta * ( i + lbnd(1) - 1 ) - r_el = sqrt ( one / ( sin(theta)**2 / initial_a**2 + & - cos(theta)**2 / initial_c**2 ) ) - xp(1) = r_el * sin(theta) + translate_x - xp(2) = translate_y - xp(3) = r_el * cos(theta) + translate_z + r_el = sqrt ( one / ( sin(theta)**2 / initial_a(l)**2 + & + cos(theta)**2 / initial_c(l)**2 ) ) + xp(1) = r_el * sin(theta) + translate_x(l) + xp(2) = translate_y(l) + xp(3) = r_el * cos(theta) + translate_z(l) xpt = matmul ( txyz, xp ) xg(i,l) = xpt(1) yg(i,l) = xpt(2) @@ -170,8 +170,8 @@ subroutine EHFinder_Init_F(CCTK_ARGUMENTS) ! if an ovaloid of Cassini is requested... if ( CCTK_EQUALS( initial_f, 'cassini' ) ) then - f(:,:,:,l) = (x**2+y**2+z**2)**2 + cas_a**4 - & - 2*cas_a**2*(x**2 - (y**2+z**2)) - cas_b**4 + f(:,:,:,l) = (x**2+y**2+z**2)**2 + cas_a(l)**4 - & + 2*cas_a(l)**2*(x**2 - (y**2+z**2)) - cas_b(l)**4 end if end do |