Skip to content

Commit

Permalink
Merge pull request #18 from ORNL-Fusion/sparse
Browse files Browse the repository at this point in the history 
Sparse
  • Loading branch information
@cianciosa
cianciosa authored Feb 7, 2024
2 parents 4158787 + 405c340 commit bf7d7c8
Show file tree
Hide file tree
Showing 27 changed files with 1,119 additions and 697 deletions.
42 changes: 21 additions & 21 deletions Sources/check3d.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -323,27 +323,27 @@ SUBROUTINE CheckForces(xc, gc)
DO ntype = 1, ndims
DO n = -ntor, ntor
DO m = 0, mpol
icol = icol+1

xc(icol,js) = eps
CALL funct_island
wplus = wtotal

xc(icol,js) = -eps
CALL funct_island
wmins = wtotal
force = -(wplus-wmins)/(2*eps)

CALL ASSERT(l_getwmhd,'L_GETWMHD = FALSE!')

xc(icol,js) = 0

IF (lWrite) THEN
nt1 = ntype
IF (ntype .GT. 3) nt1 = ntype-3
WRITE (3000+iam,100) m, n, ntype, force, force_array(nt1), &
gc1(icol,js)
END IF
icol = icol+1
xc(icol,js) = eps
CALL funct_island
wplus = wtotal

xc(icol,js) = -eps
CALL funct_island
wmins = wtotal
force = -(wplus-wmins)/(2*eps)
CALL ASSERT(l_getwmhd,'L_GETWMHD = FALSE!')

xc(icol,js) = 0

IF (lWrite) THEN
nt1 = ntype
IF (ntype .GT. 3) nt1 = ntype-3
WRITE (3000+iam,100) m, n, ntype, force, force_array(nt1), &
gc1(icol,js)
END IF

END DO
END DO
Expand Down
55 changes: 19 additions & 36 deletions Sources/diagnostics_mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -94,25 +94,6 @@ SUBROUTINE divb(ns_min, ns_max)
f_cos, nsmin, nsmax)
END IF

! DO m = 0, mpol
! DO n = -ntor, ntor
! divbmnsf(m,n,nl:nh) = ohs*(jbsupsmnsh(m,n,nl+1:nh+1) - &
! jbsupsmnsh(m,n,nl:nh)) &
! - m*(jbsupumnch(m,n,nl+1:nh+1) + &
! jbsupumnch(m,n,nl:nh))*0.5 &
! - n*nfp*(jbsupvmnch(m,n,nl+1:nh+1) + &
! jbsupvmnch(m,n,nl:nh))*0.5
! IF (lasym) THEN
! divbmncf(m,n,nl:nh) = ohs*(jbsupsmnch(m,n,nl+1:nh+1) - &
! jbsupsmnch(m,n,nl:nh)) &
! - m*(jbsupumnsh(m,n,nl+1:nh+1) + &
! jbsupumnsh(m,n,nl:nh))*0.5 &
! + n*nfp*(jbsupvmnsh(m,n,nl+1:nh+1) + &
! jbsupvmnsh(m,n,nl:nh))*0.5
! END IF
! END DO
! END DO

tnorm = hs_i*SUM(jbsupumnch(:,:,nsmin:nsmax)**2 &
+ jbsupvmnch(:,:,nsmin:nsmax)**2 &
+ jbsupsmnsh(:,:,nsmin:nsmax)**2)
Expand Down Expand Up @@ -144,7 +125,7 @@ SUBROUTINE divb(ns_min, ns_max)
!> @param[in] ns_max Upper radial index.
!-------------------------------------------------------------------------------
FUNCTION divdb(ns_min, ns_max)
USE island_params, ONLY: ohs=>ohs_i
USE island_params, ONLY: ohs=>ohs_i, fourier_context
USE hessian, ONLY: gather_array

! Declare Arguments
Expand All @@ -156,6 +137,7 @@ FUNCTION divdb(ns_min, ns_max)
INTEGER :: istat
INTEGER :: m
INTEGER :: n
INTEGER :: n_mode
INTEGER :: nl
INTEGER :: nh
REAL (dp) :: tnorm
Expand Down Expand Up @@ -184,21 +166,22 @@ FUNCTION divdb(ns_min, ns_max)
divbmnsf = 0
divbmncf = 0

DO m = 0, mpol
DO n = -ntor, ntor
divbmnsf(m,n,nl:nh) = -m*(djbsupumnch(m,n,nl+1:nh+1) + &
djbsupumnch(m,n,nl:nh))*0.5 &
- n*nfp*(djbsupvmnch(m,n,nl+1:nh+1) + &
djbsupvmnch(m,n,nl:nh))*0.5 &
+ ohs*(djbsupsmnsh(m,n,nl+1:nh+1) - &
djbsupsmnsh(m,n,nl:nh))
DO n = -ntor, ntor
n_mode = fourier_context%tor_modes(n)*nfp
DO m = 0, mpol
divbmnsf(m,n,nl:nh) = -m*(djbsupumnch(m,n,nl+1:nh+1) + &
& djbsupumnch(m,n,nl:nh))*0.5 &
& - n_mode*(djbsupvmnch(m,n,nl+1:nh+1) + &
& djbsupvmnch(m,n,nl:nh))*0.5 &
& + ohs*(djbsupsmnsh(m,n,nl+1:nh+1) - &
& djbsupsmnsh(m,n,nl:nh))
IF (lasym) THEN
divbmncf(m,n,nl:nh) = m*(djbsupumnsh(m,n,nl+1:nh+1) + &
djbsupumnsh(m,n,nl:nh))*0.5 &
+ n*nfp*(djbsupvmnsh(m,n,nl+1:nh+1) + &
djbsupvmnsh(m,n,nl:nh))*0.5 &
+ ohs*(djbsupsmnch(m,n,nl+1:nh+1) - &
djbsupsmnch(m,n,nl:nh))
divbmncf(m,n,nl:nh) = m*(djbsupumnsh(m,n,nl+1:nh+1) + &
& djbsupumnsh(m,n,nl:nh))*0.5 &
& + n_mode*(djbsupvmnsh(m,n,nl+1:nh+1) + &
& djbsupvmnsh(m,n,nl:nh))*0.5 &
& + ohs*(djbsupsmnch(m,n,nl+1:nh+1) - &
& djbsupsmnch(m,n,nl:nh))
END IF
END DO
END DO
Expand Down Expand Up @@ -428,7 +411,6 @@ SUBROUTINE SIESTA_PROFILES (iunit, arr_value, label, fsq_total)
DO n = -ntor, ntor
DO m = 0, mpol
m_array(m,n) = m
n_array(m,n) = n
END DO
END DO

Expand All @@ -437,7 +419,8 @@ SUBROUTINE SIESTA_PROFILES (iunit, arr_value, label, fsq_total)
WRITE (p_format,1000) mnmax
WRITE (iunit,p_format) ' ',' ','MPOL--->', m_array
WRITE (p_format,1000) mnmax
WRITE (iunit,p_format) 'RADIUS','RMS','NTOR--->', n_array
WRITE (iunit,p_format) 'RADIUS','RMS','NTOR--->', &
& fourier_context%tor_modes

WRITE (p_format,1002) mnmax
DO js = 2, ns
Expand Down
6 changes: 4 additions & 2 deletions Sources/evolution.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -393,6 +393,7 @@ SUBROUTINE evolve
REAL (dp), DIMENSION(1), PARAMETER :: levscan = (/ 0.25_dp /)
REAL (dp), PARAMETER :: levmarq_min = 1.0E-10_dp
REAL (dp), PARAMETER :: fsq_max = 1.E-6_dp
REAL (dp), PARAMETER :: ftol_min = 1.0E-20_dp

! Start of executable code.
l_init_state = .true.
Expand All @@ -409,8 +410,9 @@ SUBROUTINE evolve
! Determine levenberg parameter by a linear scale on how far away from the
! desired force tolarance we are.
IF (fsq_total1 .lt. fsq_max) THEN
f1 = (levmarq_param0 - levmarq_min)/(fsq_max - ftol)
f2 = (fsq_max*levmarq_min - ftol*levmarq_param0)/(fsq_max - ftol)
f1 = (levmarq_param0 - levmarq_min)/(fsq_max - MAX(ftol, ftol_min))
f2 = (fsq_max*levmarq_min - MAX(ftol, ftol_min)*levmarq_param0) &
& / (fsq_max - MAX(ftol, ftol_min))
levmarq_param = f1*fsq_total1 + f2

IF (fsq_last .le. fsq_total1) THEN
Expand Down
Loading

0 comments on commit bf7d7c8

Please sign in to comment.