Skip to content

Commit

Permalink
Feature/fortran api init functions (#165)
Browse files Browse the repository at this point in the history
* add Fortran API for MDPI initial functions

* add a fortran unit test for Fortran API for util/function/MDPI_functions
  • Loading branch information
sbrdar authored Nov 28, 2023
1 parent e864605 commit b3c4918
Show file tree
Hide file tree
Showing 7 changed files with 202 additions and 0 deletions.
15 changes: 15 additions & 0 deletions src/atlas/util/function/MDPI_functions.cc
Original file line number Diff line number Diff line change
Expand Up @@ -118,6 +118,21 @@ double MDPI_gulfstream(double lon, double lat) {
return background_func + dc * (std::max(1000. * std::sin(0.4 * (0.5 * dr + dth) + 0.007 * std::cos(50. * dth) + 0.37 * M_PI), 999.) - 999.);
}

extern "C" {
double atlas__functions__MDPI_sinusoid(double& lon, double& lat) {
return MDPI_sinusoid(lon, lat);
}
double atlas__functions__MDPI_harmonic(double& lon, double& lat) {
return MDPI_harmonic(lon, lat);
}
double atlas__functions__MDPI_vortex(double& lon, double& lat) {
return MDPI_vortex(lon, lat);
}
double atlas__functions__MDPI_gulfstream(double& lon, double& lat) {
return MDPI_gulfstream(lon, lat);
}
}

} // namespace function
} // namespace util
} // namespace atlas
7 changes: 7 additions & 0 deletions src/atlas/util/function/MDPI_functions.h
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,13 @@ double MDPI_harmonic(double lon, double lat);
double MDPI_vortex(double lon, double lat);
double MDPI_gulfstream(double lon, double lat);

extern "C" {
double atlas__functions__MDPI_sinusoid(double& lon, double& lat);
double atlas__functions__MDPI_harmonic(double& lon, double& lat);
double atlas__functions__MDPI_vortex(double& lon, double& lat);
double atlas__functions__MDPI_gulfstream(double& lon, double& lat);
}

} // namespace function

} // namespace util
Expand Down
1 change: 1 addition & 0 deletions src/atlas_f/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -205,6 +205,7 @@ ecbuild_add_library( TARGET atlas_f
SOURCES
${FORTRAN_BINDINGS}
atlas_module.F90
util/atlas_functions_module.F90
util/atlas_kinds_module.F90
util/atlas_JSON_module.F90
util/atlas_Config_module.F90
Expand Down
1 change: 1 addition & 0 deletions src/atlas_f/atlas_module.F90
Original file line number Diff line number Diff line change
Expand Up @@ -150,6 +150,7 @@ module atlas_module
& atlas_output_Gmsh
use atlas_trace_module, only : &
& atlas_Trace
use atlas_functions_module

use fckit_log_module, only: atlas_log => fckit_log

Expand Down
83 changes: 83 additions & 0 deletions src/atlas_f/util/atlas_functions_module.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,83 @@
! (C) Copyright 2013 ECMWF.
!
! This software is licensed under the terms of the Apache Licence Version 2.0
! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
! In applying this licence, ECMWF does not waive the privileges and immunities
! granted to it by virtue of its status as an intergovernmental organisation nor
! does it submit to any jurisdiction.

#include "atlas/atlas_f.h"

module atlas_functions_module

use, intrinsic :: iso_c_binding

implicit none

interface

pure function atlas__functions__MDPI_sinusoid( lon, lat ) bind(C,name="atlas__functi&
&ons__MDPI_sinusoid")
use iso_c_binding, only: c_double
real(c_double) :: atlas__functions__MDPI_sinusoid
real(c_double), intent(in) :: lon, lat
end function

pure function atlas__functions__MDPI_harmonic( lon, lat ) bind(C,name="atlas__functi&
&ons__MDPI_harmonic")
use iso_c_binding, only: c_double
real(c_double) :: atlas__functions__MDPI_harmonic
real(c_double), intent(in) :: lon, lat
end function

pure function atlas__functions__MDPI_vortex( lon, lat ) bind(C,name="atlas__function&
&s__MDPI_vortex")
use iso_c_binding, only: c_double
real(c_double) :: atlas__functions__MDPI_vortex
real(c_double), intent(in) :: lon, lat
end function

pure function atlas__functions__MDPI_gulfstream( lon, lat ) bind(C,name="atlas__func&
&tions__MDPI_gulfstream")
use iso_c_binding, only: c_double
real(c_double) :: atlas__functions__MDPI_gulfstream
real(c_double), intent(in) :: lon, lat
end function

end interface

contains

elemental function MDPI_sinusoid(lon, lat) result(val)
real(c_double), intent(in) :: lon, lat
real(c_double) :: val
val = atlas__functions__MDPI_sinusoid(lon, lat)
end function MDPI_sinusoid

! -----------------------------------------------------------------------------

elemental function MDPI_harmonic(lon, lat) result(val)
real(c_double), intent(in) :: lon, lat
real(c_double) :: val
val = atlas__functions__MDPI_harmonic(lon, lat)
end function MDPI_harmonic

! -----------------------------------------------------------------------------

elemental function MDPI_vortex(lon, lat) result(val)
real(c_double), intent(in) :: lon, lat
real(c_double) :: val
val = atlas__functions__MDPI_vortex(lon, lat)
end function MDPI_vortex

! -----------------------------------------------------------------------------

elemental function MDPI_gulfstream(lon, lat) result(val)
real(c_double), intent(in) :: lon, lat
real(c_double) :: val
val = atlas__functions__MDPI_gulfstream(lon, lat)
end function MDPI_gulfstream

! -----------------------------------------------------------------------------

end module atlas_functions_module
7 changes: 7 additions & 0 deletions src/tests/util/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,13 @@

if( HAVE_FCTEST )

add_fctest( TARGET atlas_fctest_functions
LINKER_LANGUAGE Fortran
SOURCES fctest_functions.F90
LIBS atlas_f
ENVIRONMENT ${ATLAS_TEST_ENVIRONMENT}
)

add_fctest( TARGET atlas_fctest_logging
LINKER_LANGUAGE Fortran
SOURCES fctest_logging.F90
Expand Down
88 changes: 88 additions & 0 deletions src/tests/util/fctest_functions.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@
! (C) Copyright 2013 ECMWF.
! This software is licensed under the terms of the Apache Licence Version 2.0
! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
! In applying this licence, ECMWF does not waive the privileges and immunities
! granted to it by virtue of its status as an intergovernmental organisation nor
! does it submit to any jurisdiction.

! This File contains Unit Tests for testing the
! C++ / Fortran Interfaces to the logging facilities
! @author Willem Deconinck

#include "fckit/fctest.h"

! -----------------------------------------------------------------------------

module fcta_functions_fxt
use atlas_module
use, intrinsic :: iso_c_binding
implicit none
character(len=1024) :: msg

end module fcta_functions_fxt

! -----------------------------------------------------------------------------

TESTSUITE_WITH_FIXTURE(fctest_atlas_functions,fcta_functions_fxt)

! -----------------------------------------------------------------------------
TESTSUITE_INIT
call atlas_library%initialise()
END_TESTSUITE_INIT

! -----------------------------------------------------------------------------

TESTSUITE_FINALIZE
call atlas_library%finalise()
END_TESTSUITE_FINALIZE

! -----------------------------------------------------------------------------

TEST( test_atlas_Functions )
real(c_double) :: val
val = MDPI_sinusoid(1._c_double, 1._c_double)
FCTEST_CHECK_CLOSE(val, 1.0002115216773033_c_double, 1e-12_c_double)
val = MDPI_harmonic(1._c_double, 1._c_double)
FCTEST_CHECK_CLOSE(val, 2.0000000000000000_c_double, 1e-12_c_double)
val = MDPI_vortex(1._c_double, 1._c_double)
FCTEST_CHECK_CLOSE(val, 2.7267489215500755_c_double, 1e-12_c_double)
val = MDPI_gulfstream(1._c_double, 1._c_double)
FCTEST_CHECK_CLOSE(val, 1.0002115216773033_c_double, 1e-12_c_double)
END_TEST

TEST( test_atlas_Functions_vector )
real(c_double), dimension(3) :: val, lon, lat, val_ref
lon = [ 1._c_double, 2._c_double, 3._c_double ]
lat = [ 1._c_double, 2._c_double, 3._c_double ]
val = MDPI_sinusoid(lon, lat)
val_ref = [1.0002115216773033_c_double, 1.0008458683590891_c_double, 1.0019023851484181_c_double]
FCTEST_CHECK_CLOSE(val, val_ref, 1e-12_c_double)
val = MDPI_harmonic(lon, lat)
val_ref = [2.0000000000000000_c_double, 2.0000000000000000_c_double, 2.0000000000000000_c_double]
FCTEST_CHECK_CLOSE(val, val_ref, 1e-12_c_double)
val = MDPI_vortex(lon, lat)
val_ref = [2.7267489215500755_c_double, 2.7520839004022091_c_double, 2.7755506683886928_c_double]
FCTEST_CHECK_CLOSE(val, val_ref, 1e-12_c_double)
val = MDPI_gulfstream(lon, lat)
val_ref = [1.0002115216773033_c_double, 1.0008458683590891_c_double, 1.0019023851484181_c_double]
FCTEST_CHECK_CLOSE(val, val_ref, 1e-12_c_double)
END_TEST

TEST( test_initialise_field )
type(atlas_Field) :: field_xy, field_val
real(c_double), dimension(:,:), pointer :: field_xy_v
real(c_double), dimension(:), pointer :: field_val_v
real(c_double), dimension(3) :: field_val_ref
field_xy = atlas_Field(kind=atlas_real(c_double), shape=[2,3])
field_val = atlas_Field(kind=atlas_real(c_double), shape=[3])
call field_xy%data(field_xy_v)
field_xy_v = 1._c_double
call field_val%data(field_val_v)
field_val_v = MDPI_sinusoid(field_xy_v(1,:), field_xy_v(2,:))
field_val_ref = [1.0002115216773033_c_double, 1.0002115216773033_c_double, 1.0002115216773033_c_double]
FCTEST_CHECK_CLOSE(field_val_v, field_val_ref, 1e-12_c_double)
END_TEST

! -----------------------------------------------------------------------------

END_TESTSUITE

0 comments on commit b3c4918

Please sign in to comment.