diff --git a/examples/geometry_fexample.F90 b/examples/geometry_fexample.F90 index 6c110f9..e678aab 100644 --- a/examples/geometry_fexample.F90 +++ b/examples/geometry_fexample.F90 @@ -13,21 +13,21 @@ program geometry_fexample use, intrinsic :: iso_c_binding + use Fgismo implicit none -# include "gsCInterface/gismo.ifc" character(len=80, kind=C_CHAR) :: some_file - type(c_ptr) :: g + type(t_gsgeometry) :: g integer :: nu, nv ! TODO: input option for xml-file? write(*,'(2(a,f5.1),a,i3)') 'reading XML for tensor B-spline' some_file = 'sw_tp.xml' // C_NULL_CHAR - g = gsCReadFile(some_file) + call f_gsCReadFile(some_file, g) - write(*,'(a,i3)') 'done, g.dim=', gsFunctionSet_domainDim(g) + write(*,'(a,i3)') 'done, g.dim=', f_gsFunctionSet_domainDim(g) - call gsFunctionSet_print(g) + call f_gsFunctionSet_print(g) if (.true.) then call show_basic_usage( g ) @@ -37,7 +37,7 @@ program geometry_fexample call show_recover_points( g ) endif - call gsFunctionSet_delete(g) + call f_gsfunctionset_delete(g) write(*,*) 'done.' end program geometry_fexample @@ -49,9 +49,8 @@ subroutine show_basic_usage( g ) use, intrinsic :: iso_c_binding use Fgismo implicit none -# include "gsCInterface/gismo.ifc" !--subroutine arguments - type(c_ptr) :: g + type(t_gsgeometry) :: g !--local variables integer(C_INT) :: nRows, nCols, out_rows, out_cols, irow, icol, icoor, ipar type(t_gsmatrix) :: uvm, xyzm @@ -75,8 +74,7 @@ subroutine show_basic_usage( g ) uvm = f_gsmatrix_create_rcd(nRows, nCols, uv) xyzm = f_gsmatrix_create() - call gsFunctionSet_eval_into(G, uvm%c_mat, xyzm%c_mat) - call f_gsmatrix_update_data_ptr( xyzm ) + call f_gsFunctionSet_eval_into(G, uvm, xyzm) ! call f_gsmatrix_print(xyzm) ! show output data @@ -102,9 +100,8 @@ subroutine show_recover_points( g ) use, intrinsic :: iso_c_binding use Fgismo implicit none -# include "gsCInterface/gismo.ifc" !--subroutine arguments - type(c_ptr) :: g + type(t_gsgeometry) :: g !--local variables integer(C_INT), parameter :: XDIR = 0, YDIR = 1, ZDIR = 2 integer(C_INT) :: nCols, irow, icol, out_rows, out_cols @@ -135,9 +132,7 @@ subroutine show_recover_points( g ) uvm = f_gsmatrix_create() eps = 1d-6 - call gsGeometry_recoverPoints(G, uvm%c_mat, xyzm%c_mat, ZDIR, eps) - call f_gsmatrix_update_data_ptr( uvm ) - call f_gsmatrix_update_data_ptr( xyzm ) + call f_gsGeometry_recoverPoints(G, uvm, xyzm, ZDIR, eps) ! print output data @@ -162,4 +157,3 @@ subroutine show_recover_points( g ) end subroutine show_recover_points !----------------------------------------------------------------------------------------------------------- - diff --git a/src/Fgismo.F90 b/src/Fgismo.F90 index 5e6fd7c..db8d7f6 100644 --- a/src/Fgismo.F90 +++ b/src/Fgismo.F90 @@ -4,20 +4,13 @@ module Fgismo implicit none private - interface -# include "gsCMatrix.ifc" -# include "gsCMatrixInt.ifc" -# include "gsCVector.ifc" -# include "gsCVectorInt.ifc" -# include "gsCKnotVector.ifc" -# include "gsCFunctionSet.ifc" -# include "gsCMultiPatch.ifc" -# include "gsCBasis.ifc" -# include "gsCGeometry.ifc" -# include "gsCReadFile.ifc" - end interface +# include "gismo.ifc" public t_gsmatrix + public t_gsfunctionset + public t_gsgeometry + + ! functions from gsCMatrix.ifc: public f_gsmatrix_create public f_gsmatrix_create_rcd public f_gsmatrix_print @@ -27,6 +20,18 @@ module Fgismo public f_gsmatrix_delete public f_gsmatrix_update_data_ptr + ! functions from gsCFunctionSet.ifc: + public f_gsfunctionset_delete + public f_gsfunctionset_print + public f_gsfunctionset_domaindim + public f_gsfunctionset_eval_into + + ! functions from gsCReadFile.ifc: + public f_gscreadfile + + ! functions from gsCGeometry.ifc: + public f_gsgeometry_recoverpoints + !------------------------------------------------------------------------------------------------------------ type :: t_gsmatrix @@ -34,44 +39,55 @@ module Fgismo real(C_DOUBLE), dimension(:,:), pointer :: data => NULL() ! link to array in C/C++ gsMatrix end type t_gsmatrix +!------------------------------------------------------------------------------------------------------------ + + type :: t_gsfunctionset + type(C_PTR) :: c_fs ! C/C++ gsFunctionSet object + end type t_gsfunctionset + +!------------------------------------------------------------------------------------------------------------ + + type, extends(t_gsfunctionset) :: t_gsgeometry ! C/C++ gsGeometry object + end type t_gsgeometry + !------------------------------------------------------------------------------------------------------------ contains +!------------------------------------------------------------------------------------------------------------ +! wrap functions from gsCMatrix.ifc: !------------------------------------------------------------------------------------------------------------ -function f_gsmatrix_create() +function f_gsmatrix_create() result(f_mat) !--purpose: create empty gsmatrix object #ifdef _WIN32 !dec$ attributes dllexport :: f_gsmatrix_create #endif implicit none -!--function result type: - type(t_gsmatrix) :: f_gsmatrix_create +!--function result: + type(t_gsmatrix) :: f_mat - f_gsmatrix_create%c_mat = gsMatrix_create() - call f_gsmatrix_update_data_ptr( f_gsmatrix_create ) + f_mat%c_mat = gsMatrix_create() + call f_gsmatrix_update_data_ptr( f_mat ) end function f_gsmatrix_create !------------------------------------------------------------------------------------------------------------ -function f_gsmatrix_create_rcd(nrows, ncols, data) +function f_gsmatrix_create_rcd(nrows, ncols, data) result(f_mat) !--purpose: create gsmatrix object from input data #ifdef _WIN32 !dec$ attributes dllexport :: f_gsmatrix_create_rcd #endif implicit none -!--function result type: - type(t_gsmatrix) :: f_gsmatrix_create_rcd +!--function result: + type(t_gsmatrix) :: f_mat !--function arguments: integer(C_INT) :: nrows, ncols real(C_DOUBLE), dimension(:,:) :: data - associate( f_mat => f_gsmatrix_create_rcd ) f_mat%c_mat = gsMatrix_create_rcd(nrows, ncols, data) call f_gsmatrix_update_data_ptr( f_mat ) - end associate end function f_gsmatrix_create_rcd @@ -92,34 +108,34 @@ end subroutine f_gsmatrix_print !------------------------------------------------------------------------------------------------------------ -function f_gsmatrix_rows(f_mat) +function f_gsmatrix_rows(f_mat) result(nrows) !--purpose: get number of rows from gsmatrix object #ifdef _WIN32 !dec$ attributes dllexport :: f_gsmatrix_rows #endif implicit none -!--function return value: - integer(C_INT) :: f_gsmatrix_rows +!--function result: + integer(C_INT) :: nrows !--function arguments: type(t_gsmatrix) :: f_mat - f_gsmatrix_rows = gsMatrix_rows(f_mat%c_mat) + nrows = gsMatrix_rows(f_mat%c_mat) end function f_gsmatrix_rows !------------------------------------------------------------------------------------------------------------ -function f_gsmatrix_cols(f_mat) +function f_gsmatrix_cols(f_mat) result(ncols) !--purpose: get number of columns from gsmatrix object #ifdef _WIN32 !dec$ attributes dllexport :: f_gsmatrix_cols #endif implicit none -!--function return value: - integer(C_INT) :: f_gsmatrix_cols +!--function result: + integer(C_INT) :: ncols !--function arguments: type(t_gsmatrix) :: f_mat - f_gsmatrix_cols = gsMatrix_cols(f_mat%c_mat) + ncols = gsMatrix_cols(f_mat%c_mat) end function f_gsmatrix_cols !------------------------------------------------------------------------------------------------------------ @@ -146,19 +162,19 @@ end subroutine f_gsmatrix_update_data_ptr !------------------------------------------------------------------------------------------------------------ -function f_gsmatrix_data(f_mat) +function f_gsmatrix_data(f_mat) result(f_data_ptr) !--purpose: get pointer to data array of a gsmatrix object #ifdef _WIN32 !dec$ attributes dllexport :: f_gsmatrix_data #endif implicit none -!--function return value: - real(C_DOUBLE), dimension(:,:), pointer :: f_gsmatrix_data +!--function result: + real(C_DOUBLE), dimension(:,:), pointer :: f_data_ptr !--function arguments: type(t_gsmatrix) :: f_mat call f_gsmatrix_update_data_ptr(f_mat) - f_gsmatrix_data => f_mat%data + f_data_ptr => f_mat%data end function f_gsmatrix_data !------------------------------------------------------------------------------------------------------------ @@ -175,6 +191,114 @@ subroutine f_gsmatrix_delete(f_mat) nullify(f_mat%data) end subroutine f_gsmatrix_delete +!------------------------------------------------------------------------------------------------------------ +! wrap functions of gsFunctionSet.ifc: +!------------------------------------------------------------------------------------------------------------ + +subroutine f_gsfunctionset_delete(f_fs) +!--purpose: destroy a gsfunctionset object +#ifdef _WIN32 +!dec$ attributes dllexport :: f_gsfunctionset_delete +#endif +!--subroutine arguments: + class(t_gsfunctionset) :: f_fs + + call gsFunctionSet_delete(f_fs%c_fs) +end subroutine f_gsfunctionset_delete + +!------------------------------------------------------------------------------------------------------------ + +subroutine f_gsfunctionset_print(f_fs) +!--purpose: print contents of a gsfunctionset object +#ifdef _WIN32 +!dec$ attributes dllexport :: f_gsfunctionset_print +#endif +!--subroutine arguments: + class(t_gsfunctionset) :: f_fs + + call gsfunctionset_print(f_fs%c_fs) +end subroutine f_gsfunctionset_print + +!------------------------------------------------------------------------------------------------------------ + +function f_gsfunctionset_domaindim(f_fs) result(domdim) +!--purpose: get domain dimension of a gsfunctionset object +#ifdef _WIN32 +!dec$ attributes dllexport :: f_gsfunctionset_domaindim +#endif +!--function result: + integer(C_INT) :: domdim +!--function arguments: + class(t_gsfunctionset) :: f_fs + + domdim = gsfunctionset_domaindim(f_fs%c_fs) +end function f_gsfunctionset_domaindim + +!------------------------------------------------------------------------------------------------------------ +! GISMO_EXPORT int gsFunctionSet_targetDim(gsCFunctionSet * fs); + +!------------------------------------------------------------------------------------------------------------ +! GISMO_EXPORT gsCMatrix* gsFunctionSet_support(gsCFunctionSet * fs); + +!------------------------------------------------------------------------------------------------------------ + +subroutine f_gsfunctionset_eval_into(f_fs, f_uv, f_result) +!--purpose: evaluate gsfunctionset object at parameter values uv into result matrix +#ifdef _WIN32 +!dec$ attributes dllexport :: f_gsfunctionset_eval_into +#endif +!--subroutine arguments: + class(t_gsfunctionset) :: f_fs + type(t_gsmatrix) :: f_uv, f_result + + call gsfunctionset_eval_into(f_fs%c_fs, f_uv%c_mat, f_result%c_mat ) + call f_gsmatrix_update_data_ptr( f_result ) +end subroutine f_gsfunctionset_eval_into + +!------------------------------------------------------------------------------------------------------------ +! GISMO_EXPORT void gsFunctionSet_deriv_into(gsCFunctionSet * fs, +! gsCMatrix * u, +! gsCMatrix * result); + +!------------------------------------------------------------------------------------------------------------ + +!------------------------------------------------------------------------------------------------------------ +! wrap functions of gsCReadFile.ifc: +!------------------------------------------------------------------------------------------------------------ + +subroutine f_gscreadfile(filename, f_fs) +!--purpose: read a gsgeometry object from a file +#ifdef _WIN32 +!dec$ attributes dllexport :: f_gscreadfile +#endif +!--subroutine arguments: + character(len=1,kind=C_CHAR) :: filename(*) + class(t_gsfunctionset) :: f_fs + + f_fs%c_fs = gsCReadFile(filename) +end subroutine f_gscreadfile + +!------------------------------------------------------------------------------------------------------------ +! wrap functions of gsCGeometry.ifc: +!------------------------------------------------------------------------------------------------------------ + +subroutine f_gsgeometry_recoverpoints(f_geom, f_uv, f_xyz, k, accuracy) +!--purpose: invert coordinates (x,y), (x,z), or (y,z) to corresponding (u,v) and fill in missing coordinate +#ifdef _WIN32 +!dec$ attributes dllexport :: f_gsgeometry_recoverpoints +#endif +!--subroutine arguments: + type(t_gsgeometry) :: f_geom + type(t_gsmatrix) :: f_uv, f_xyz + integer(C_INT) :: k + real(C_DOUBLE) :: accuracy + + call gsGeometry_recoverPoints(f_geom%c_fs, f_uv%c_mat, f_xyz%c_mat, k, accuracy) + call f_gsmatrix_update_data_ptr( f_uv ) + call f_gsmatrix_update_data_ptr( f_xyz ) + +end subroutine f_gsgeometry_recoverpoints + !------------------------------------------------------------------------------------------------------------ end module Fgismo