Skip to content

Commit bcc0ef9

Browse files
committed
[fitpack] implement parallel comm routines
1 parent 77259fc commit bcc0ef9

File tree

2 files changed

+297
-0
lines changed

2 files changed

+297
-0
lines changed

src/fitpack_core.f90

Lines changed: 198 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,7 @@
1919
! **************************************************************************************************
2020
module fitpack_core
2121
use iso_c_binding, only: c_double,c_int32_t,c_bool
22+
use iso_fortran_env, only: real64, int32
2223
implicit none
2324
private
2425

@@ -27,6 +28,10 @@ module fitpack_core
2728
integer, parameter, public :: FP_SIZE = c_int32_t
2829
integer, parameter, public :: FP_FLAG = c_int32_t
2930
integer, parameter, public :: FP_BOOL = c_bool
31+
integer, parameter, public :: FP_COMM = c_double
32+
33+
!> Marker for unallocated arrays in communication buffers
34+
integer(FP_SIZE), parameter :: FP_NOT_ALLOC = -99999999_FP_SIZE
3035

3136
! Curve fitting routines
3237
public :: curfit ! * General curve fitting
@@ -79,6 +84,12 @@ module fitpack_core
7984
public :: fitpack_error_handling
8085
public :: get_smoothing
8186

87+
! Parallel communication utilities for standalone use
88+
public :: FP_RCOMMS_PER_BITS
89+
public :: FP_COMM_SIZE
90+
public :: FP_COMM_PACK
91+
public :: FP_COMM_EXPAND
92+
8293
! Spline behavior for points not in the support
8394
integer(FP_FLAG), parameter, public :: OUTSIDE_EXTRAPOLATE = 0 ! extrapolated from the end spans
8495
integer(FP_FLAG), parameter, public :: OUTSIDE_ZERO = 1 ! spline evaluates to zero
@@ -167,6 +178,27 @@ end function fitpack_polar_boundary
167178
module procedure swap_size
168179
end interface fitpack_swap
169180

181+
!> Communication size for allocatable 1D arrays
182+
interface FP_COMM_SIZE
183+
module procedure FP_REAL_COMM_SIZE_1D
184+
module procedure FP_REAL_COMM_SIZE_2D
185+
module procedure FP_SIZE_COMM_SIZE_1D
186+
end interface FP_COMM_SIZE
187+
188+
!> Pack allocatable 1D arrays into communication buffer
189+
interface FP_COMM_PACK
190+
module procedure FP_REAL_COMM_PACK_1D
191+
module procedure FP_REAL_COMM_PACK_2D
192+
module procedure FP_SIZE_COMM_PACK_1D
193+
end interface FP_COMM_PACK
194+
195+
!> Expand communication buffer into allocatable 1D arrays
196+
interface FP_COMM_EXPAND
197+
module procedure FP_REAL_COMM_EXPAND_1D
198+
module procedure FP_REAL_COMM_EXPAND_2D
199+
module procedure FP_SIZE_COMM_EXPAND_1D
200+
end interface FP_COMM_EXPAND
201+
170202
contains
171203

172204
! Flow control: on output flag present, return it;
@@ -18726,4 +18758,170 @@ elemental logical(FP_BOOL) function not_equal(a,b)
1872618758
not_equal = .not.equal(a,b)
1872718759
end function not_equal
1872818760

18761+
! ==============================================================================================
18762+
! PARALLEL COMMUNICATION UTILITIES
18763+
! ==============================================================================================
18764+
18765+
!> Number of FP_COMM elements required to store nbits of data
18766+
elemental integer(FP_SIZE) function FP_RCOMMS_PER_BITS(nbits)
18767+
integer(FP_SIZE), intent(in) :: nbits
18768+
integer, parameter :: COMM_BITS = storage_size(0.0_FP_COMM)
18769+
FP_RCOMMS_PER_BITS = nbits / COMM_BITS + merge(IONE, IZERO, mod(nbits, COMM_BITS) > 0)
18770+
end function FP_RCOMMS_PER_BITS
18771+
18772+
!> Calculate storage size for 1D real(FP_REAL) allocatable array
18773+
!> Header: 1 FP_COMM storing bounds as 2 int32 (or NOT_ALLOC marker)
18774+
!> Data: raw transfer of array contents
18775+
pure integer(FP_SIZE) function FP_REAL_COMM_SIZE_1D(array)
18776+
real(FP_REAL), allocatable, intent(in) :: array(:)
18777+
integer(FP_SIZE) :: n
18778+
18779+
FP_REAL_COMM_SIZE_1D = rank(array) ! Header
18780+
if (allocated(array)) then
18781+
n = size(array, kind=FP_SIZE)
18782+
FP_REAL_COMM_SIZE_1D = FP_REAL_COMM_SIZE_1D + FP_RCOMMS_PER_BITS(n * storage_size(array))
18783+
end if
18784+
end function FP_REAL_COMM_SIZE_1D
18785+
18786+
!> Calculate storage size for 2D real(FP_REAL) allocatable array
18787+
pure integer(FP_SIZE) function FP_REAL_COMM_SIZE_2D(array)
18788+
real(FP_REAL), allocatable, intent(in) :: array(:,:)
18789+
integer(FP_SIZE) :: n
18790+
18791+
FP_REAL_COMM_SIZE_2D = rank(array) ! Header (2 dimensions)
18792+
if (allocated(array)) then
18793+
n = size(array, kind=FP_SIZE)
18794+
FP_REAL_COMM_SIZE_2D = FP_REAL_COMM_SIZE_2D + FP_RCOMMS_PER_BITS(n * storage_size(array))
18795+
end if
18796+
end function FP_REAL_COMM_SIZE_2D
18797+
18798+
!> Calculate storage size for 1D integer(FP_SIZE) allocatable array
18799+
pure integer(FP_SIZE) function FP_SIZE_COMM_SIZE_1D(array)
18800+
integer(FP_SIZE), allocatable, intent(in) :: array(:)
18801+
integer(FP_SIZE) :: n
18802+
18803+
FP_SIZE_COMM_SIZE_1D = rank(array) ! Header
18804+
if (allocated(array)) then
18805+
n = size(array, kind=FP_SIZE)
18806+
FP_SIZE_COMM_SIZE_1D = FP_SIZE_COMM_SIZE_1D + FP_RCOMMS_PER_BITS(n * storage_size(array))
18807+
end if
18808+
end function FP_SIZE_COMM_SIZE_1D
18809+
18810+
!> Pack 1D real(FP_REAL) allocatable array into communication buffer
18811+
pure subroutine FP_REAL_COMM_PACK_1D(array, buffer)
18812+
real(FP_REAL), allocatable, intent(in) :: array(:)
18813+
real(FP_COMM), intent(out) :: buffer(:)
18814+
18815+
integer(FP_SIZE) :: bnd(2), ndoubles
18816+
integer(FP_SIZE), parameter :: header = 1
18817+
18818+
if (allocated(array)) then
18819+
bnd = [lbound(array, 1, FP_SIZE), ubound(array, 1, FP_SIZE)]
18820+
else
18821+
bnd = FP_NOT_ALLOC
18822+
end if
18823+
18824+
buffer(1:header) = transfer(bnd, buffer(1:header), int(header))
18825+
18826+
if (all(bnd /= FP_NOT_ALLOC)) then
18827+
ndoubles = FP_RCOMMS_PER_BITS(size(array, kind=FP_SIZE) * storage_size(array))
18828+
buffer(header+1:header+ndoubles) = transfer(array, buffer(header+1:header+ndoubles), int(ndoubles))
18829+
end if
18830+
end subroutine FP_REAL_COMM_PACK_1D
18831+
18832+
!> Pack 2D real(FP_REAL) allocatable array into communication buffer
18833+
pure subroutine FP_REAL_COMM_PACK_2D(array, buffer)
18834+
real(FP_REAL), allocatable, intent(in) :: array(:,:)
18835+
real(FP_COMM), intent(out) :: buffer(:)
18836+
18837+
integer(FP_SIZE) :: bnd(2, 2), ndoubles, d
18838+
integer(FP_SIZE), parameter :: header = 2
18839+
18840+
if (allocated(array)) then
18841+
forall (d = 1:2) bnd(:, d) = [lbound(array, d, FP_SIZE), ubound(array, d, FP_SIZE)]
18842+
else
18843+
bnd = FP_NOT_ALLOC
18844+
end if
18845+
18846+
buffer(1:header) = transfer(bnd, buffer(1:header), int(header))
18847+
18848+
if (all(bnd /= FP_NOT_ALLOC)) then
18849+
ndoubles = FP_RCOMMS_PER_BITS(size(array, kind=FP_SIZE) * storage_size(array))
18850+
buffer(header+1:header+ndoubles) = transfer(array, buffer(header+1:header+ndoubles), int(ndoubles))
18851+
end if
18852+
end subroutine FP_REAL_COMM_PACK_2D
18853+
18854+
!> Pack 1D integer(FP_SIZE) allocatable array into communication buffer
18855+
pure subroutine FP_SIZE_COMM_PACK_1D(array, buffer)
18856+
integer(FP_SIZE), allocatable, intent(in) :: array(:)
18857+
real(FP_COMM), intent(out) :: buffer(:)
18858+
18859+
integer(FP_SIZE) :: bnd(2), ndoubles
18860+
integer(FP_SIZE), parameter :: header = 1
18861+
18862+
if (allocated(array)) then
18863+
bnd = [lbound(array, 1, FP_SIZE), ubound(array, 1, FP_SIZE)]
18864+
else
18865+
bnd = FP_NOT_ALLOC
18866+
end if
18867+
18868+
buffer(1:header) = transfer(bnd, buffer(1:header), int(header))
18869+
18870+
if (all(bnd /= FP_NOT_ALLOC)) then
18871+
ndoubles = FP_RCOMMS_PER_BITS(size(array, kind=FP_SIZE) * storage_size(array))
18872+
buffer(header+1:header+ndoubles) = transfer(array, buffer(header+1:header+ndoubles), int(ndoubles))
18873+
end if
18874+
end subroutine FP_SIZE_COMM_PACK_1D
18875+
18876+
!> Expand communication buffer into 1D real(FP_REAL) allocatable array
18877+
pure subroutine FP_REAL_COMM_EXPAND_1D(array, buffer)
18878+
real(FP_REAL), allocatable, intent(out) :: array(:)
18879+
real(FP_COMM), intent(in) :: buffer(:)
18880+
18881+
integer(FP_SIZE) :: bnd(2), n
18882+
integer(FP_SIZE), parameter :: header = 1
18883+
18884+
bnd = transfer(buffer(:header), bnd)
18885+
18886+
if (all(bnd /= FP_NOT_ALLOC)) then
18887+
allocate(array(bnd(1):bnd(2)))
18888+
n = FP_RCOMMS_PER_BITS(size(array, kind=FP_SIZE) * storage_size(array))
18889+
array = transfer(buffer(header+1:header+n), array)
18890+
end if
18891+
end subroutine FP_REAL_COMM_EXPAND_1D
18892+
18893+
!> Expand communication buffer into 2D real(FP_REAL) allocatable array
18894+
pure subroutine FP_REAL_COMM_EXPAND_2D(array, buffer)
18895+
real(FP_REAL), allocatable, intent(out) :: array(:,:)
18896+
real(FP_COMM), intent(in) :: buffer(:)
18897+
18898+
integer(FP_SIZE) :: bnd(2, 2), n
18899+
integer(FP_SIZE), parameter :: header = 2
18900+
18901+
bnd = reshape(transfer(buffer(:header), bnd), shape(bnd))
18902+
18903+
if (all(bnd /= FP_NOT_ALLOC)) then
18904+
allocate(array(bnd(1,1):bnd(2,1), bnd(1,2):bnd(2,2)))
18905+
n = FP_RCOMMS_PER_BITS(size(array, kind=FP_SIZE) * storage_size(array))
18906+
array = reshape(transfer(buffer(header+1:header+n), array, size(array)), shape(array))
18907+
end if
18908+
end subroutine FP_REAL_COMM_EXPAND_2D
18909+
18910+
!> Expand communication buffer into 1D integer(FP_SIZE) allocatable array
18911+
pure subroutine FP_SIZE_COMM_EXPAND_1D(array, buffer)
18912+
integer(FP_SIZE), allocatable, intent(out) :: array(:)
18913+
real(FP_REAL), intent(in) :: buffer(:)
18914+
18915+
integer(FP_SIZE) :: bnd(2), ndoubles
18916+
integer(FP_SIZE), parameter :: header = 1
18917+
18918+
bnd = transfer(buffer(:header), bnd)
18919+
18920+
if (all(bnd /= FP_NOT_ALLOC)) then
18921+
allocate(array(bnd(1):bnd(2)))
18922+
ndoubles = FP_RCOMMS_PER_BITS(size(array, kind=FP_SIZE) * storage_size(array))
18923+
array = transfer(buffer(header+1:header+ndoubles), array)
18924+
end if
18925+
end subroutine FP_SIZE_COMM_EXPAND_1D
18926+
1872918927
end module fitpack_core

src/fitpack_curves.f90

Lines changed: 99 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -112,6 +112,11 @@ module fitpack_curves
112112
!> Properties: MSE
113113
procedure, non_overridable :: mse => curve_error
114114

115+
!> Parallel communication interface (size/pack/expand)
116+
procedure :: comm_size => curve_comm_size
117+
procedure :: comm_pack => curve_comm_pack
118+
procedure :: comm_expand => curve_comm_expand
119+
115120
end type fitpack_curve
116121

117122
!> Derived type describing a periodic curve. No changes are made to the storage,
@@ -618,4 +623,98 @@ function zeros(this,ierr)
618623

619624
end function zeros
620625

626+
! =================================================================================================
627+
! PARALLEL COMMUNICATION (size/pack/expand)
628+
! =================================================================================================
629+
630+
!> Return communication buffer size (number of FP_REAL elements)
631+
!> This counts storage for all curve data needed to reconstruct the spline
632+
pure integer(FP_SIZE) function curve_comm_size(this)
633+
class(fitpack_curve), intent(in) :: this
634+
635+
! Scalar integers: m, order, knots, bc, iopt, nest, lwrk (7 values)
636+
! Scalar reals: xleft, xright, smoothing, fp (4 values)
637+
! Use 11 FP_REAL slots for 11 scalar values
638+
curve_comm_size = 11 &
639+
+ FP_COMM_SIZE(this%x) &
640+
+ FP_COMM_SIZE(this%y) &
641+
+ FP_COMM_SIZE(this%sp) &
642+
+ FP_COMM_SIZE(this%w) &
643+
+ FP_COMM_SIZE(this%t) &
644+
+ FP_COMM_SIZE(this%c) &
645+
+ FP_COMM_SIZE(this%wrk) &
646+
+ FP_COMM_SIZE(this%iwrk) &
647+
+ FP_COMM_SIZE(this%wrk_fou)
648+
649+
end function curve_comm_size
650+
651+
!> Pack curve data into communication buffer
652+
pure subroutine curve_comm_pack(this, buffer)
653+
class(fitpack_curve), intent(in) :: this
654+
real(FP_REAL), intent(out) :: buffer(:)
655+
656+
integer(FP_SIZE) :: pos, n
657+
658+
! Pack scalar integers as reals
659+
buffer(1) = real(this%m, FP_REAL)
660+
buffer(2) = real(this%order, FP_REAL)
661+
buffer(3) = real(this%knots, FP_REAL)
662+
buffer(4) = real(this%bc, FP_REAL)
663+
buffer(5) = real(this%iopt, FP_REAL)
664+
buffer(6) = real(this%nest, FP_REAL)
665+
buffer(7) = real(this%lwrk, FP_REAL)
666+
buffer(8) = this%xleft
667+
buffer(9) = this%xright
668+
buffer(10) = this%smoothing
669+
buffer(11) = this%fp
670+
pos = 12
671+
672+
! Pack 1D real arrays using FP_COMM_PACK
673+
call FP_COMM_PACK(this%x, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%x)
674+
call FP_COMM_PACK(this%y, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%y)
675+
call FP_COMM_PACK(this%sp, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%sp)
676+
call FP_COMM_PACK(this%w, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%w)
677+
call FP_COMM_PACK(this%t, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%t)
678+
call FP_COMM_PACK(this%c, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%c)
679+
call FP_COMM_PACK(this%wrk, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%wrk)
680+
call FP_COMM_PACK(this%iwrk, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%iwrk)
681+
call FP_COMM_PACK(this%wrk_fou, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%wrk_fou)
682+
683+
end subroutine curve_comm_pack
684+
685+
!> Expand curve data from communication buffer
686+
pure subroutine curve_comm_expand(this, buffer)
687+
class(fitpack_curve), intent(inout) :: this
688+
real(FP_REAL), intent(in) :: buffer(:)
689+
690+
integer(FP_SIZE) :: pos
691+
692+
! Expand scalar integers
693+
this%m = nint(buffer(1), FP_SIZE)
694+
this%order = nint(buffer(2), FP_SIZE)
695+
this%knots = nint(buffer(3), FP_SIZE)
696+
this%bc = nint(buffer(4), FP_FLAG)
697+
this%iopt = nint(buffer(5), FP_SIZE)
698+
this%nest = nint(buffer(6), FP_SIZE)
699+
this%lwrk = nint(buffer(7), FP_SIZE)
700+
this%xleft = buffer(8)
701+
this%xright = buffer(9)
702+
this%smoothing = buffer(10)
703+
this%fp = buffer(11)
704+
pos = 12
705+
706+
! Expand arrays using FP_COMM_EXPAND
707+
! After expand, FP_COMM_SIZE returns size based on new allocation
708+
call FP_COMM_EXPAND(this%x, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%x)
709+
call FP_COMM_EXPAND(this%y, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%y)
710+
call FP_COMM_EXPAND(this%sp, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%sp)
711+
call FP_COMM_EXPAND(this%w, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%w)
712+
call FP_COMM_EXPAND(this%t, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%t)
713+
call FP_COMM_EXPAND(this%c, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%c)
714+
call FP_COMM_EXPAND(this%wrk, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%wrk)
715+
call FP_COMM_EXPAND(this%iwrk, buffer(pos:)); pos = pos + FP_COMM_SIZE(this%iwrk)
716+
call FP_COMM_EXPAND(this%wrk_fou, buffer(pos:))
717+
718+
end subroutine curve_comm_expand
719+
621720
end module fitpack_curves

0 commit comments

Comments
 (0)