Skip to content

Commit

Permalink
240125.170942.HKT C interface: add comments to algorithm_c.f90, add d…
Browse files Browse the repository at this point in the history
…eallocation
  • Loading branch information
zaikunzhang committed Jan 25, 2024
1 parent 5ba7c90 commit 9880e58
Show file tree
Hide file tree
Showing 5 changed files with 189 additions and 115 deletions.
63 changes: 40 additions & 23 deletions c/bobyqa_c.f90
Original file line number Diff line number Diff line change
Expand Up @@ -27,8 +27,7 @@ subroutine bobyqa_c(cobj_ptr, data_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, &
type(C_FUNPTR), intent(IN), value :: cobj_ptr
type(C_PTR), intent(in), value :: data_ptr
integer(C_INT), intent(in), value :: n
! We cannot use assumed-shape arrays for C interoperability
real(C_DOUBLE), intent(inout) :: x(n)
real(C_DOUBLE), intent(inout) :: x(n) ! We cannot use assumed-shape arrays for C interoperability
real(C_DOUBLE), intent(out) :: f
type(C_PTR), intent(in), value :: xl
type(C_PTR), intent(in), value :: xu
Expand All @@ -45,34 +44,45 @@ subroutine bobyqa_c(cobj_ptr, data_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, &
! Local variables
integer(IK) :: info_loc
integer(IK) :: iprint_loc
integer(IK) :: nf_loc
integer(IK), allocatable :: maxfun_loc
integer(IK), allocatable :: npt_loc
integer(IK) :: nf_loc
! The initialization below to null is necessary to avoid a bug with the newer Intel compiler ifx.
! See details here: https://fortran-lang.discourse.group/t/strange-issue-with-ifx-compiler-and-assume-recursion/7013
! The bug was observed in all versions of ifx up to 2024.0.1. Once this bug is fixed we should remove the
! initialization to null because it implies the 'save' attribute, which is undesirable.
procedure(CCALLBACK), pointer :: cb_ptr => null()
procedure(COBJ), pointer :: obj_ptr => null()
real(C_DOUBLE), pointer :: xl_loc_interm(:)
real(C_DOUBLE), pointer :: xu_loc_interm(:)
real(RP) :: f_loc
real(RP), allocatable :: rhobeg_loc
real(RP), allocatable :: rhoend_loc
real(RP) :: ftarget_loc
real(RP) :: x_loc(n)
real(C_DOUBLE), pointer :: xl_loc_interm(:)
real(RP), allocatable :: rhobeg_loc
real(RP), allocatable :: rhoend_loc
real(RP), allocatable :: xl_loc(:)
real(C_DOUBLE), pointer :: xu_loc_interm(:)
real(RP), allocatable :: xu_loc(:)
! The initialization to null is necessary to avoid a bug with the newer Intel compiler ifx.
! See details here: https://fortran-lang.discourse.group/t/strange-issue-with-ifx-compiler-and-assume-recursion/7013
! The bug was observed in all versions of ifx up to 2024.0.1. Once this bug is fixed we should remove the
! initialization to null because it implies the 'save' attribute, which is undesirable.
procedure(COBJ), pointer :: obj_ptr => null()
procedure(CCALLBACK), pointer :: cb_ptr => null()

! Read the inputs and convert them to the Fortran side types

! The following inputs correspond to compulsory arguments in the Fortran code.
x_loc = real(x, kind(x_loc))
if (C_ASSOCIATED(xl)) then
call C_F_POINTER(xl, xl_loc_interm, shape=[n])
call c_f_procpointer(cobj_ptr, obj_ptr)

! The following inputs correspond to optional arguments in the Fortran code.
! Since C does not support optional arguments, we use NaN to represent an absent real scalar, 0 to
! represent an absent integer scalar (all integer arguments are expected positive), and an
! unassociated pointer to represent an absent array. In case of NaN, 0, and unassociated pointers,
! the allocatable variables such as RHOBEG_LOC will be left uninitialized and hence unallocated, and
! then treated as an absent argument when passed to the Fortran code.
! See Sec. 9.7.1.3 and 15.5.2.13 of J3/24-007 (Fortran 2023 Interpretation Document).
if (c_associated(xl)) then
call c_f_pointer(xl, xl_loc_interm, shape=[n])
call safealloc(xl_loc, int(n, IK))
xl_loc = real(xl_loc_interm, kind(xl_loc))
end if
if (C_ASSOCIATED(xu)) then
call C_F_POINTER(xu, xu_loc_interm, shape=[n])
if (c_associated(xu)) then
call c_f_pointer(xu, xu_loc_interm, shape=[n])
call safealloc(xu_loc, int(n, IK))
xu_loc = real(xu_loc_interm, kind(xu_loc))
end if
Expand All @@ -90,13 +100,12 @@ subroutine bobyqa_c(cobj_ptr, data_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, &
npt_loc = int(npt, kind(npt_loc))
end if
iprint_loc = int(iprint, kind(iprint_loc))
call C_F_PROCPOINTER(cobj_ptr, obj_ptr)

! Call the Fortran code
if (C_ASSOCIATED(callback_ptr)) then
if (c_associated(callback_ptr)) then
! If a C callback function is provided, we convert it to a Fortran procedure pointer and capture
! that pointer in the closure below.
call C_F_PROCPOINTER(callback_ptr, cb_ptr)
call c_f_procpointer(callback_ptr, cb_ptr)
! We then provide the closure to the algorithm.
call bobyqa(calfun, x_loc, f_loc, xl=xl_loc, xu=xu_loc, nf=nf_loc, rhobeg=rhobeg_loc, rhoend=rhoend_loc, &
& ftarget=ftarget_loc, maxfun=maxfun_loc, npt=npt_loc, iprint=iprint_loc, callback_fcn=callback_fcn, info=info_loc)
Expand All @@ -111,6 +120,14 @@ subroutine bobyqa_c(cobj_ptr, data_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, &
nf = int(nf_loc, kind(nf))
info = int(info_loc, kind(info))

! Deallocate variables not needed any more. Indeed, automatic allocation will take place at exit.
if (allocated(npt_loc)) deallocate (npt_loc)
if (allocated(maxfun_loc)) deallocate (maxfun_loc)
if (allocated(rhoend_loc)) deallocate (rhoend_loc)
if (allocated(rhobeg_loc)) deallocate (rhobeg_loc)
if (allocated(xu_loc)) deallocate (xu_loc)
if (allocated(xu_loc)) deallocate (xl_loc)

contains

!--------------------------------------------------------------------------------------------------!
Expand Down Expand Up @@ -201,12 +218,12 @@ subroutine callback_fcn(x_sub, f_sub, nf_sub, tr, cstrv_sub, nlconstr_sub, termi
call cb_ptr(n_sub_loc, x_sub_loc, f_sub_loc, nf_sub_loc, tr_loc, cstrv_sub_loc, m_nlconstr, nlconstr_sub_loc, terminate_loc)

! Write the output
if ( present(terminate) ) then
if (present(terminate)) then
terminate = logical(terminate_loc, kind(terminate))
end if

! Deallocate resources
if (allocated(nlconstr_sub_loc)) deallocate(nlconstr_sub_loc)
! Deallocate variables not needed any more. Indeed, automatic allocation will take place at exit.
if (allocated(nlconstr_sub_loc)) deallocate (nlconstr_sub_loc)

end subroutine callback_fcn

Expand Down
81 changes: 49 additions & 32 deletions c/cobyla_c.f90
Original file line number Diff line number Diff line change
Expand Up @@ -27,8 +27,7 @@ subroutine cobyla_c(m_nlcon, cobjcon_ptr, data_ptr, n, x, f, cstrv, nlconstr, m_
type(C_FUNPTR), intent(in), value :: cobjcon_ptr
type(C_PTR), intent(in), value :: data_ptr
integer(C_INT), intent(in), value :: n
! We cannot use assumed-shape arrays for C interoperability
real(C_DOUBLE), intent(inout) :: x(n)
real(C_DOUBLE), intent(inout) :: x(n) ! We cannot use assumed-shape arrays for C interoperability
real(C_DOUBLE), intent(out) :: f
real(C_DOUBLE), intent(out) :: cstrv
real(C_DOUBLE), intent(out) :: nlconstr(m_nlcon)
Expand All @@ -55,53 +54,65 @@ subroutine cobyla_c(m_nlcon, cobjcon_ptr, data_ptr, n, x, f, cstrv, nlconstr, m_
integer(IK) :: info_loc
integer(IK) :: iprint_loc
integer(IK) :: m_nlcon_loc
integer(IK), allocatable :: maxfun_loc
integer(IK) :: nf_loc
real(RP) :: Aineq_loc(m_ineq, n)
real(RP) :: bineq_loc(m_ineq)
integer(IK), allocatable :: maxfun_loc
! The initialization below to null is necessary to avoid a bug with the newer Intel compiler ifx.
! See details here: https://fortran-lang.discourse.group/t/strange-issue-with-ifx-compiler-and-assume-recursion/7013
! The bug was observed in all versions of ifx up to 2024.0.1. Once this bug is fixed we should remove the
! initialization to null because it implies the 'save' attribute, which is undesirable.
procedure(CCALLBACK), pointer :: cb_ptr => null()
procedure(COBJCON), pointer :: objcon_ptr => null()
real(C_DOUBLE), pointer :: nlconstr0_loc_interm(:)
real(C_DOUBLE), pointer :: xl_loc_interm(:)
real(C_DOUBLE), pointer :: xu_loc_interm(:)
real(RP) :: Aeq_loc(m_eq, n)
real(RP) :: Aineq_loc(m_ineq, n)
real(RP) :: beq_loc(m_eq)
real(RP) :: bineq_loc(m_ineq)
real(RP) :: cstrv_loc
real(RP) :: nlconstr_loc(m_nlcon)
real(RP) :: f_loc
real(RP), allocatable :: rhobeg_loc
real(RP), allocatable :: rhoend_loc
real(RP) :: ftarget_loc
real(RP) :: nlconstr_loc(m_nlcon)
real(RP) :: x_loc(n)
real(C_DOUBLE), pointer :: xl_loc_interm(:)
real(RP), allocatable :: xl_loc(:)
real(C_DOUBLE), pointer :: xu_loc_interm(:)
real(RP), allocatable :: xu_loc(:)
real(RP), allocatable :: f0_loc
real(C_DOUBLE), pointer :: nlconstr0_loc_interm(:)
real(RP), allocatable :: nlconstr0_loc(:)
! The initialization to null is necessary to avoid a bug with the newer Intel compiler ifx.
! See details here: https://fortran-lang.discourse.group/t/strange-issue-with-ifx-compiler-and-assume-recursion/7013
! The bug was observed in all versions of ifx up to 2024.0.1. Once this bug is fixed we should remove the
! initialization to null because it implies the 'save' attribute, which is undesirable.
procedure(COBJCON), pointer :: objcon_ptr => null()
procedure(CCALLBACK), pointer :: cb_ptr => null()
real(RP), allocatable :: rhobeg_loc
real(RP), allocatable :: rhoend_loc
real(RP), allocatable :: xl_loc(:)
real(RP), allocatable :: xu_loc(:)

! Read the inputs and convert them to the Fortran side types
! Note that `transpose` is needed when reading 2D arrays, since they are stored in the row-major
! order in c but column-major in Fortran.

! The following inputs correspond to compulsory arguments in the Fortran code.
x_loc = real(x, kind(x_loc))
m_nlcon_loc = int(m_nlcon, kind(m_nlcon_loc))
call c_f_procpointer(cobjcon_ptr, objcon_ptr)

! The following inputs correspond to optional arguments in the Fortran code.
! Since C does not support optional arguments, we use NaN to represent an absent real scalar, 0 to
! represent an absent integer scalar (all integer arguments are expected positive), and an
! unassociated pointer to represent an absent array. In case of NaN, 0, and unassociated pointers,
! the allocatable variables such as RHOBEG_LOC will be left uninitialized and hence unallocated, and
! then treated as an absent argument when passed to the Fortran code.
! See Sec. 9.7.1.3 and 15.5.2.13 of J3/24-007 (Fortran 2023 Interpretation Document).
Aineq_loc = real(transpose(Aineq), kind(Aineq_loc))
bineq_loc = real(bineq, kind(bineq_loc))
Aeq_loc = real(transpose(Aeq), kind(Aeq_loc))
beq_loc = real(beq, kind(beq_loc))
if (C_ASSOCIATED(xl)) then
call C_F_POINTER(xl, xl_loc_interm, shape=[n])
if (c_associated(xl)) then
call c_f_pointer(xl, xl_loc_interm, shape=[n])
call safealloc(xl_loc, int(n, IK))
xl_loc = real(xl_loc_interm, kind(xl_loc))
end if
if (C_ASSOCIATED(xu)) then
call C_F_POINTER(xu, xu_loc_interm, shape=[n])
if (c_associated(xu)) then
call c_f_pointer(xu, xu_loc_interm, shape=[n])
call safealloc(xu_loc, int(n, IK))
xu_loc = real(xu_loc_interm, kind(xu_loc))
end if
if (C_ASSOCIATED(nlconstr0)) then
call C_F_POINTER(nlconstr0, nlconstr0_loc_interm, shape=[m_nlcon])
if (c_associated(nlconstr0)) then
call c_f_pointer(nlconstr0, nlconstr0_loc_interm, shape=[m_nlcon])
call safealloc(nlconstr0_loc, int(m_nlcon, IK))
nlconstr0_loc = real(nlconstr0_loc_interm, kind(nlconstr0_loc))
! We assume that if nlconstr0 was provided, that the f0 provided is valid
Expand All @@ -118,14 +129,12 @@ subroutine cobyla_c(m_nlcon, cobjcon_ptr, data_ptr, n, x, f, cstrv, nlconstr, m_
maxfun_loc = int(maxfun, kind(maxfun_loc))
end if
iprint_loc = int(iprint, kind(iprint_loc))
m_nlcon_loc = int(m_nlcon, kind(m_nlcon_loc))
call C_F_PROCPOINTER(cobjcon_ptr, objcon_ptr)

! Call the Fortran code
if (C_ASSOCIATED(callback_ptr)) then
if (c_associated(callback_ptr)) then
! If a C callback function is provided, we convert it to a Fortran procedure pointer and capture
! that pointer in the closure below.
call C_F_PROCPOINTER(callback_ptr, cb_ptr)
call c_f_procpointer(callback_ptr, cb_ptr)
! We then provide the closure to the algorithm.
call cobyla(calcfc, m_nlcon_loc, x_loc, f_loc, cstrv=cstrv_loc, nlconstr=nlconstr_loc, Aineq=Aineq_loc, &
& bineq=bineq_loc, Aeq=Aeq_loc, beq=beq_loc, xl=xl_loc, xu=xu_loc, f0=f0_loc, nlconstr0=nlconstr0_loc, &
Expand All @@ -146,6 +155,14 @@ subroutine cobyla_c(m_nlcon, cobjcon_ptr, data_ptr, n, x, f, cstrv, nlconstr, m_
info = int(info_loc, kind(info))
nlconstr = real(nlconstr_loc, kind(nlconstr))

! Deallocate variables not needed any more. Indeed, automatic allocation will take place at exit.
if (allocated(npt_loc)) deallocate (npt_loc)
if (allocated(maxfun_loc)) deallocate (maxfun_loc)
if (allocated(rhoend_loc)) deallocate (rhoend_loc)
if (allocated(rhobeg_loc)) deallocate (rhobeg_loc)
if (allocated(xu_loc)) deallocate (xu_loc)
if (allocated(xu_loc)) deallocate (xl_loc)

contains

!--------------------------------------------------------------------------------------------------!
Expand Down Expand Up @@ -239,12 +256,12 @@ subroutine callback_fcn(x_sub, f_sub, nf_sub, tr, cstrv_sub, nlconstr_sub, termi
call cb_ptr(n_sub_loc, x_sub_loc, f_sub_loc, nf_sub_loc, tr_loc, cstrv_sub_loc, m_nlconstr, nlconstr_sub_loc, terminate_loc)

! Write the output
if ( present(terminate) ) then
if (present(terminate)) then
terminate = logical(terminate_loc, kind(terminate))
end if

! Deallocate resources
if (allocated(nlconstr_sub_loc)) deallocate(nlconstr_sub_loc)
! Deallocate variables not needed any more. Indeed, automatic allocation will take place at exit.
if (allocated(nlconstr_sub_loc)) deallocate (nlconstr_sub_loc)

end subroutine callback_fcn

Expand Down
Loading

1 comment on commit 9880e58

@github-actions
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@check-spelling-bot Report

🔴 Please review

See the 📜action log or 📝 job summary for details.

Unrecognized words (1)

procpointer

To accept these unrecognized words as correct, you could run the following commands

... in a clone of the [email protected]:libprima/prima.git repository
on the main branch (ℹ️ how do I use this?):

curl -s -S -L 'https://raw.githubusercontent.com/check-spelling/check-spelling/main/apply.pl' |
perl - 'https://github.com/libprima/prima/actions/runs/7652241197/attempts/1'
If the flagged items are 🤯 false positives

If items relate to a ...

  • binary file (or some other file you wouldn't want to check at all).

    Please add a file path to the excludes.txt file matching the containing file.

    File paths are Perl 5 Regular Expressions - you can test yours before committing to verify it will match your files.

    ^ refers to the file's path from the root of the repository, so ^README\.md$ would exclude README.md (on whichever branch you're using).

  • well-formed pattern.

    If you can write a pattern that would match it,
    try adding it to the patterns.txt file.

    Patterns are Perl 5 Regular Expressions - you can test yours before committing to verify it will match your lines.

    Note that patterns can't match multiline strings.

Please sign in to comment.