From b939f96d64d17545fb0d8d499261a2675f73a03d Mon Sep 17 00:00:00 2001
From: cianciosa <m4c@ornl.gov>
Date: Fri, 4 Oct 2024 14:26:06 -0400
Subject: [PATCH] If free boundary restart. There was a unnecessary
normalization factor being applied.
---
Sources/grid_extension.f90 | 20 ++++++--------------
Sources/metrics.f90 | 2 +-
Sources/pchelms.f90 | 4 ++--
3 files changed, 9 insertions(+), 17 deletions(-)
diff --git a/Sources/grid_extension.f90 b/Sources/grid_extension.f90
index 0566bc4..0919983 100644
--- a/Sources/grid_extension.f90
+++ b/Sources/grid_extension.f90
@@ -143,8 +143,6 @@ SUBROUTINE grid_extender(wout_file, itype, istat)
CLASS (fourier_class), POINTER :: four => null()
! Start of executable code.
- IF (.NOT.l_vessel) RETURN ! FIXME: This check should be outsize the subroutine
-
CALL init_extender
! chfit = 'quad'
@@ -193,8 +191,8 @@ SUBROUTINE grid_extender(wout_file, itype, istat)
! Compute extended "s" coordinate for s>1, assuming Vp(s) = Vp(s=1)*s
! where Vp(s=1) is obtained from extrapolating VMEC values. Results is
! V(s) = Vp*(s^2-1)/2 + V0 => s_edge = SQRT(2(V1-V0)/Vp + 1)
- vol0 = Get_Volume(zero) !rho = 0: s=1 vmec
- vol1 = Get_Volume(one) !rho = 1: vessel
+ vol0 = get_volume(zero) !rho = 0: s=1 vmec
+ vol1 = get_volume(one) !rho = 1: vessel
vp1 = GetVP1()
s_edge = SQRT(1 + 2*(vol1 - vol0)/vp1)
@@ -255,7 +253,7 @@ SUBROUTINE grid_extender(wout_file, itype, istat)
END IF
END DO
- nse = nse-1 !Not to double count s=1 surface
+ nse = nse - 1 !Not to double count s=1 surface
ALLOCATE(r_full(ntheta,nzeta,nse), z_full(ntheta,nzeta,nse), &
ru_full(ntheta,nzeta,nse), zu_full(ntheta,nzeta,nse), &
rv_full(ntheta,nzeta,nse), zv_full(ntheta,nzeta,nse))
@@ -324,12 +322,6 @@ SUBROUTINE grid_extender(wout_file, itype, istat)
f_cos)
CALL fourier_context%tomnsp(z1_i(:,:,nsp1:ns), zmns_i(:,:,nsp1:ns), &
f_sin)
- DO js = nsp1, ns
- rmnc_i(:,:,js) = rmnc_i(:,:,js) &
- / fourier_context%orthonorm(0:mpol,-ntor:ntor)
- zmns_i(:,:,js) = zmns_i(:,:,js) &
- / fourier_context%orthonorm(0:mpol,-ntor:ntor)
- END DO
! PUT PROFILES ONTO EXTENDED MESH
tmp = 0
@@ -403,9 +395,9 @@ SUBROUTINE grid_extender(wout_file, itype, istat)
END IF
#if defined(NOSKS)
- DO u=1, ntheta
+ DO js=1, ns
DO v=1, nzeta
- DO js=1, ns
+ DO u=1, ntheta
WRITE(1600,*) r1_i(u,v,js), z1_i(u,v,js)
CALL FLUSH(1600)
END DO
@@ -434,7 +426,7 @@ SUBROUTINE grid_extender(wout_file, itype, istat)
DEALLOCATE (r12, rs12, zs12, ru12, zu12)
- ming = MINVAL(sqrtg_full(:,:,ns:nse) )
+ ming = MINVAL(sqrtg_full(:,:,ns:nse))
maxg = MAXVAL(sqrtg_full(:,:,ns:nse))
PRINT *,' IAM: ', iam, 'MING: ', ming,' MAXG: ',maxg
CALL ASSERT(ming*maxg.GT.zero,' Jacobian changed sign in grid_extender')
diff --git a/Sources/metrics.f90 b/Sources/metrics.f90
index f051955..259678d 100644
--- a/Sources/metrics.f90
+++ b/Sources/metrics.f90
@@ -375,7 +375,7 @@ SUBROUTINE half_mesh_metrics(r1_i, ru_i, rv_i, z1_i, zu_i, zv_i)
mintest = MINVAL(sqrtg(:,2:))
maxtest = MAXVAL(sqrtg(:,2:))
-
+
IF (mintest*maxtest .LE. zero) THEN
imin = MINLOC(sqrtg(:,2:))
imax = MAXLOC(sqrtg(:,2:))
diff --git a/Sources/pchelms.f90 b/Sources/pchelms.f90
index e59ac53..2d22ae2 100644
--- a/Sources/pchelms.f90
+++ b/Sources/pchelms.f90
@@ -90,7 +90,7 @@ SUBROUTINE run_pchelms
LOGICAL :: lcolscale_save
! local parameters
- REAL (dp), PARAMETER :: levmarq_param_inital = 1.E-6_dp
+! REAL (dp), PARAMETER :: levmarq_param_inital = 1.E-6_dp
! Start of executable code
siesta_curtor = 0.0
@@ -109,7 +109,7 @@ SUBROUTINE run_pchelms
CALL second0(ton)
! Add initial diffusive term to preconditioner instead
- levmarq_param = levmarq_param_inital
+ !levmarq_param = levmarq_param_inital
! Note initHess must have already been called prior to this point.
CALL Compute_Hessian_Blocks(Compute_Forces_Lin, .FALSE.)