Dear Jed, dear Barry, Thank you for the fast answers ! If I have any success I will make sure to make a pull request to provide some version of ex11 in Fortran. Regarding the stubs, I admit I started looking in that direction to add the missing wrappers but I am not sure I fully understand the process yet. For each C function, I gotta provide a Fortran interface in a .h90 file as well as a C function that has a Fortran-like prototype and calls the C function - right ? However there are a few things I could not find / understand yet. For instance, it appears that for C functions that have character string arguments take an extra argument in their Fortran-like-prototype-wrapper, namely the length of the string. Is that passed automatically ? I couldn’t find where it could come from ... Another thing is for functions like PetscFVView. I guess the wrapping is less straightforward because I tried a quick something and it segfault’ed. I couldnt find the wrapper for DMView although there is such a routine in Fortran too. Could you please detail how to wrap such functions ? Thank you very much again, Thibault Bridel-Bertomeu Le sam. 12 déc. 2020 à 21:48, Barry Smith <[email protected]> a écrit :

PETSc Fortran interfaces are a combination of automatically generated and manually generated.

For any C PETSc function if the manual page begins with /*@ it generates the Fortran interface automatically (make allfortranstubs). If it begins /*@C then either the Fortran interface is done manually or is missing.

C functions that have character string arguments or function arguments (or a few other special cases) need to be manually provided. The automatically generated stubs go in the directory ftn-auto while manually generated ones go in the directory fin-custom.

Perhaps you could first generate a list of "missing" Fortran stubs and then for each stub determine why it is missing and if it can be provided. Some are likely easy to provide but a few (involving function arguments) will be more involved. Once you have all the stubs available translating ex11.c becomes straightforward.

Barry

On Dec 12, 2020, at 9:30 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

Is there somewhere a version of the TS tutorial ex11.c in Fortran ? I am looking into building in F90 (let's say that it is an unavoidable constraint) an unstructured 3D solver of the Euler equations using the "new" features of PETSc - mostly DMPlex & PetscFV - but I think there are some interfaces missing and I find it hard to find workarounds in Fortran. I would be grateful if anyone could please give me some pointers ...

Thank you very much in advance,

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

--

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

Good morning all, Thank you Barry for your answer. I started adding some interfaces (PetscFVSetComponentName, PetscFVView & PetscFVSetType so far) and I have to say I think it is working quite well, but the prototypes of those functions are still quite "simple". I am stuck at how to implement the wrappers for PetscDSSetRiemannSolver and PetscDSSetContext though, especially on how to pass a function as an argument to PetscDSSetRiemannSolver ... Are there any similar functions that may already have their wrappers ? Thank you very much, Thibault Le sam. 12 déc. 2020 à 23:28, Barry Smith <[email protected]> a écrit :

On Dec 12, 2020, at 2:59 PM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear Jed, dear Barry,

Thank you for the fast answers !

If I have any success I will make sure to make a pull request to provide some version of ex11 in Fortran.

Regarding the stubs, I admit I started looking in that direction to add the missing wrappers but I am not sure I fully understand the process yet. For each C function, I gotta provide a Fortran interface in a .h90 file as well as a C function that has a Fortran-like prototype and calls the C function - right ?

Yes,

However there are a few things I could not find / understand yet. For instance, it appears that for C functions that have character string arguments take an extra argument in their Fortran-like-prototype-wrapper, namely the length of the string. Is that passed automatically ? I couldn’t find where it could come from ...

This secret argument is put automatically by the Fortran compiler.

Another thing is for functions like PetscFVView. I guess the wrapping is less straightforward because I tried a quick something and it segfault’ed. I couldnt find the wrapper for DMView although there is such a routine in Fortran too. Could you please detail how to wrap such functions ?

PETSC_EXTERN void dmview_(DM *da,PetscViewer *vin,PetscErrorCode *ierr) { PetscViewer v; PetscPatchDefaultViewers_Fortran(vin,v); *ierr = DMView(*da,v); }

dm/interface/ftn-custom/zdmf.c

Thank you very much again,

Thibault Bridel-Bertomeu

Le sam. 12 déc. 2020 à 21:48, Barry Smith <[email protected]> a écrit :

...

PETSc Fortran interfaces are a combination of automatically generated and manually generated.

For any C PETSc function if the manual page begins with /*@ it generates the Fortran interface automatically (make allfortranstubs). If it begins /*@C then either the Fortran interface is done manually or is missing.

C functions that have character string arguments or function arguments (or a few other special cases) need to be manually provided. The automatically generated stubs go in the directory ftn-auto while manually generated ones go in the directory fin-custom.

Perhaps you could first generate a list of "missing" Fortran stubs and then for each stub determine why it is missing and if it can be provided. Some are likely easy to provide but a few (involving function arguments) will be more involved. Once you have all the stubs available translating ex11.c becomes straightforward.

Barry

On Dec 12, 2020, at 9:30 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

Is there somewhere a version of the TS tutorial ex11.c in Fortran ? I am looking into building in F90 (let's say that it is an unavoidable constraint) an unstructured 3D solver of the Euler equations using the "new" features of PETSc - mostly DMPlex & PetscFV - but I think there are some interfaces missing and I find it hard to find workarounds in Fortran. I would be grateful if anyone could please give me some pointers ...

Thank you very much in advance,

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

--

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

Thank you Mark for your answer. I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ? On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ? Thank you ! Thibault Le dim. 13 déc. 2020 à 15:17, Mark Adams <[email protected]> a écrit :

I don't think function pointers to PETSc (DS and DM) methods are not going to work in Fortran:

ierr = PetscDSAddBoundary(prob, DM_BC_NATURAL_RIEMANN, "inflow", "Face Sets", 0, 0, NULL, (void (*)(void)) PhysicsBoundary_Advect_Inflow, NULL, ALEN(inflowids), inflowids, phys);CHKERRQ(ierr);

You could write a funcs.c file that you call from your fortran code, like, call setBC1(prob,...,ierr)

and put PhysicsBoundary_Advect_Inflow and setBC1 in funcs.c, for instance.

On Sun, Dec 13, 2020 at 5:32 AM Thibault Bridel-Bertomeu < [email protected]> wrote:

...

Good morning all,

Thank you Barry for your answer. I started adding some interfaces (PetscFVSetComponentName, PetscFVView & PetscFVSetType so far) and I have to say I think it is working quite well, but the prototypes of those functions are still quite "simple". I am stuck at how to implement the wrappers for PetscDSSetRiemannSolver and PetscDSSetContext though, especially on how to pass a function as an argument to PetscDSSetRiemannSolver ... Are there any similar functions that may already have their wrappers ?

Thank you very much,

Thibault

Le sam. 12 déc. 2020 à 23:28, Barry Smith <[email protected]> a écrit :

...

On Dec 12, 2020, at 2:59 PM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear Jed, dear Barry,

Thank you for the fast answers !

If I have any success I will make sure to make a pull request to provide some version of ex11 in Fortran.

Regarding the stubs, I admit I started looking in that direction to add the missing wrappers but I am not sure I fully understand the process yet. For each C function, I gotta provide a Fortran interface in a .h90 file as well as a C function that has a Fortran-like prototype and calls the C function - right ?

Yes,

However there are a few things I could not find / understand yet. For instance, it appears that for C functions that have character string arguments take an extra argument in their Fortran-like-prototype-wrapper, namely the length of the string. Is that passed automatically ? I couldn’t find where it could come from ...

This secret argument is put automatically by the Fortran compiler.

Another thing is for functions like PetscFVView. I guess the wrapping is less straightforward because I tried a quick something and it segfault’ed. I couldnt find the wrapper for DMView although there is such a routine in Fortran too. Could you please detail how to wrap such functions ?

PETSC_EXTERN void dmview_(DM *da,PetscViewer *vin,PetscErrorCode *ierr) { PetscViewer v; PetscPatchDefaultViewers_Fortran(vin,v); *ierr = DMView(*da,v); }

dm/interface/ftn-custom/zdmf.c

Thank you very much again,

Thibault Bridel-Bertomeu

Le sam. 12 déc. 2020 à 21:48, Barry Smith <[email protected]> a écrit :

...

PETSc Fortran interfaces are a combination of automatically generated and manually generated.

For any C PETSc function if the manual page begins with /*@ it generates the Fortran interface automatically (make allfortranstubs). If it begins /*@C then either the Fortran interface is done manually or is missing.

C functions that have character string arguments or function arguments (or a few other special cases) need to be manually provided. The automatically generated stubs go in the directory ftn-auto while manually generated ones go in the directory fin-custom.

Perhaps you could first generate a list of "missing" Fortran stubs and then for each stub determine why it is missing and if it can be provided. Some are likely easy to provide but a few (involving function arguments) will be more involved. Once you have all the stubs available translating ex11.c becomes straightforward.

Barry

On Dec 12, 2020, at 9:30 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

Is there somewhere a version of the TS tutorial ex11.c in Fortran ? I am looking into building in F90 (let's say that it is an unavoidable constraint) an unstructured 3D solver of the Euler equations using the "new" features of PETSc - mostly DMPlex & PetscFV - but I think there are some interfaces missing and I find it hard to find workarounds in Fortran. I would be grateful if anyone could please give me some pointers ...

Thank you very much in advance,

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

--

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

Hello everyone, So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already. Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread. Stay safe, Thibault Bridel-Bertomeu Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec *,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy); } else { *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx); *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

Hello everyone, Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ? Thibault Bridel-Bertomeu Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

This is great. If you make a branch off of the PETSc git repository with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec *,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy); } else { *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx); *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

Hello Barry, I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ? Thibault Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

Thibault,

A subdirectory for the example is fine; we have other examples that use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git repository with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec *,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy); } else { *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx); *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry ! The MR will appear in a little while ! Thibault Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected]> a écrit :

Good question. There is a trick to limit the amount of work you need to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected]:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

...

Thibault,

A subdirectory for the example is fine; we have other examples that use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git repository with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec *,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy); } else { *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx); *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

PetscObjectUseFortranCallback((PetscDS)ctx,

*ierr = PetscObjectSetFortranCallback((PetscObject)*prob

It looks like the problem is that these user provided functions do not take a PetscDS directly as an argument so the Fortran callback information cannot be obtained from them. The manual page for PetscDSAddBoundary() says - ctx - An optional user context for bcFunc but then when it lists the calling sequence for bcFunc it does not list the ctx as an argument, so either the manual page or code is wrong. It looks like you make the ctx be the PetscDS prob argument when you call PetscDSAddBoundary In principle this sounds like it might work. I think you need to track through the debugger to see if the ctx passed to ourbocofunc() is actually the PetscDS prob variable and if not why it is not. Barry

On Dec 22, 2020, at 5:49 AM, Thibault Bridel-Bertomeu <[email protected]> wrote:

Dear all,

I have hit two snags while implementing the missing wrappers necessary to transcribe ex11 to Fortran.

First is about the PetscDSAddBoundary wrapper, that I have done so :

static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, ctx); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, *prob); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); }

but when I call it in the program, with adequate routines, I obtain the following error :

[0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: Corrupt argument: https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind <https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind> [0]PETSC ERROR: Fortran callback not set on this object [0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html <https://www.mcs.anl.gov/petsc/documentation/faq.html> for trouble shooting. [0]PETSC ERROR: Petsc Development GIT revision: v3.14.2-297-gf36a7edeb8 GIT Date: 2020-12-18 04:42:53 +0000 [0]PETSC ERROR: ../../../bin/eulerian3D on a named macbook-pro-de-thibault.home by tbridel Sun Dec 20 15:05:15 2020 [0]PETSC ERROR: Configure options --with-clean=0 --prefix=/Users/tbridel/Documents/1-CODES/04-PETSC/build --with-make-np=2 --with-windows-graphics=0 --with-debugging=0 --download-fblaslapack --download-mpich-shared=0 --with-x=0 --with-pthread=0 --with-valgrind=0 --PETSC_ARCH=macosx --with-fc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpifort --with-cc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpicc --with-cxx=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpic++ --with-openmp=0 --download-hypre=yes --download-sowing=yes --download-metis=yes --download-parmetis=yes --download-triangle=yes --download-tetgen=yes --download-ctetgen=yes --download-p4est=yes --download-zlib=yes --download-c2html=yes --download-eigen=yes --download-pragmatic=yes --with-hdf5-dir=/usr/local/Cellar/hdf5/1.10.5_1 --with-cmake-dir=/usr/local/Cellar/cmake/3.15.3 [0]PETSC ERROR: #1 PetscObjectGetFortranCallback() line 258 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/sys/objects/inherit.c [0]PETSC ERROR: #2 ourbocofunc() line 141 in /Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c [0]PETSC ERROR: #3 DMPlexInsertBoundaryValuesRiemann() line 989 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c [0]PETSC ERROR: #4 DMPlexInsertBoundaryValues_Plex() line 1052 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c [0]PETSC ERROR: #5 DMPlexInsertBoundaryValues() line 1142 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c [0]PETSC ERROR: #6 DMPlexComputeResidual_Internal() line 4524 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c [0]PETSC ERROR: #7 DMPlexTSComputeRHSFunctionFVM() line 74 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmplexts.c [0]PETSC ERROR: #8 ourdmtsrhsfunc() line 186 in /Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c [0]PETSC ERROR: #9 TSComputeRHSFunction_DMLocal() line 105 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmlocalts.c [0]PETSC ERROR: #10 TSComputeRHSFunction() line 653 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c [0]PETSC ERROR: #11 TSSSPStep_RK_3() line 120 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c [0]PETSC ERROR: #12 TSStep_SSP() line 208 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c [0]PETSC ERROR: #13 TSStep() line 3757 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c [0]PETSC ERROR: #14 TSSolve() line 4154 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c [0]PETSC ERROR: #15 User provided function() line 0 in User file

Second is about the DMProjectFunction wrapper, that I have done so :

static PetscErrorCode ourdmprojfunc(PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) { PetscObjectUseFortranCallback((DM)ctx, dmprojfunc, (PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*), (&dim,&time,x,&Nf,u,_ctx,&ierr)) } PETSC_EXTERN void dmprojectfunction_(DM *dm, PetscReal *time, void (*func)(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*), void *ctx, InsertMode *mode, Vec X, PetscErrorCode *ierr) { PetscErrorCode (*funcarr[1]) (PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_CLASS, &dmprojfunc, (PetscVoidFunction)func, ctx); funcarr[0] = ourdmprojfunc; *ierr = DMProjectFunction(*dm, *time, funcarr, &ctx, *mode, X); }

This time there is no error because I cannot reach this point in the program, but I am not sure anyways how to write this wrapper, especially because of the double pointers that DMProjectFunction takes as arguments.

Does anyone have any idea what could be going wrong with those two wrappers ?

Thank you very much in advance !!

Thibault

Le ven. 18 déc. 2020 à 11:02, Thibault Bridel-Bertomeu <[email protected] <mailto:[email protected]>> a écrit : Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry !

The MR will appear in a little while !

Thibault

Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected] <mailto:[email protected]>> a écrit :

Good question. There is a trick to limit the amount of work you need to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected] <mailto:[email protected]>:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

...

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu <[email protected] <mailto:[email protected]>> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected] <mailto:[email protected]>> a écrit :

Thibault,

A subdirectory for the example is fine; we have other examples that use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

...

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu <[email protected] <mailto:[email protected]>> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected] <mailto:[email protected]>> a écrit :

This is great. If you make a branch off of the PETSc git repository with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/ <https://docs.petsc.org/en/latest/developers/integration/>

Barry

...

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu <[email protected] <mailto:[email protected]>> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected] <mailto:[email protected]>> a écrit : Thibault Bridel-Bertomeu <[email protected] <mailto:[email protected]>> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) { PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec *,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy); } else { *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx); *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } } <wrapper_petsc.h90><wrapper_petsc.c>

Thibault Bridel-Bertomeu <[email protected]> writes:

Good morning everyone,

Thank you Barry for the answer, it works now !

I am facing (yet) another situation: the TSAdaptRegister function. In the MR on gitlab, Jed mentioned that sometimes, when function pointers are not stored in PETSc objects, one can use stack memory to pass that pointer from fortran to C.

The issue with stack memory is that when it returns, that memory is invalid. You can't use it in this instance. I think you're going to have problems implementing a TSAdaptCreate_XYZ in Fortran (because the body of that function will need to access private struct members; see below). I would implement what you need in C and you can call out to Fortran if you want from inside TSAdaptChoose_YourMethod(). PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt adapt) { TSAdapt_DSP *dsp; PetscErrorCode ierr; PetscFunctionBegin; ierr = PetscNewLog(adapt,&dsp);CHKERRQ(ierr); adapt->reject_safety = 1.0; /* unused */ adapt->data = (void*)dsp; adapt->ops->choose = TSAdaptChoose_DSP; adapt->ops->setfromoptions = TSAdaptSetFromOptions_DSP; adapt->ops->destroy = TSAdaptDestroy_DSP; adapt->ops->view = TSAdaptView_DSP; ierr = PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetFilter_C",TSAdaptDSPSetFilter_DSP);CHKERRQ(ierr); ierr = PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetPID_C",TSAdaptDSPSetPID_DSP);CHKERRQ(ierr); ierr = TSAdaptDSPSetFilter_DSP(adapt,"PI42");CHKERRQ(ierr); ierr = TSAdaptRestart_DSP(adapt);CHKERRQ(ierr); PetscFunctionReturn(0); }

Can anyone develop that idea ? Because for TSAdaptRegister, i guess the wrapper would start like :

PETSC_EXTERN void tsadaptregister_(char *sname, void (*func)(TSAdapt*,PetscErrorCode*), PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T snamelen)

but then the C TSAdaptRegister function takes a PetscErrorCode (*func)(TSAdapt) function pointer as argument ... I cannot use any FORTRAN_CALLBACK here since I do not have any object to hook it to, and I could not find a similar situation among the pre-existing wrappers. Does anyone have an idea on how to proceed ?

Thanks !!

Thibault

Le mar. 22 déc. 2020 à 21:20, Barry Smith <[email protected]> a écrit :

...

PetscObjectUseFortranCallback((PetscDS)ctx,

*ierr = PetscObjectSetFortranCallback((PetscObject)*prob

It looks like the problem is that these user provided functions do not take a PetscDS directly as an argument so the Fortran callback information cannot be obtained from them.

The manual page for PetscDSAddBoundary() says

- ctx - An optional user context for bcFunc

but then when it lists the calling sequence for bcFunc it does not list the ctx as an argument, so either the manual page or code is wrong.

It looks like you make the ctx be the PetscDS prob argument when you call PetscDSAddBoundary

In principle this sounds like it might work. I think you need to track through the debugger to see if the ctx passed to ourbocofunc() is actually the PetscDS prob variable and if not why it is not.

Barry

On Dec 22, 2020, at 5:49 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

I have hit two snags while implementing the missing wrappers necessary to transcribe ex11 to Fortran.

First is about the PetscDSAddBoundary wrapper, that I have done so :

static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, ctx); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, *prob); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); }

but when I call it in the program, with adequate routines, I obtain the following error :

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------[0]PETSC ERROR: Corrupt argument: https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC ERROR: Fortran callback not set on this object[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.[0]PETSC ERROR: Petsc Development GIT revision: v3.14.2-297-gf36a7edeb8 GIT Date: 2020-12-18 04:42:53 +0000[0]PETSC ERROR: ../../../bin/eulerian3D on a named macbook-pro-de-thibault.home by tbridel Sun Dec 20 15:05:15 2020[0]PETSC ERROR: Configure options --with-clean=0 --prefix=/Users/tbridel/Documents/1-CODES/04-PETSC/build --with-make-np=2 --with-windows-graphics=0 --with-debugging=0 --download-fblaslapack --download-mpich-shared=0 --with-x=0 --with-pthread=0 --with-valgrind=0 --PETSC_ARCH=macosx --with-fc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpifort --with-cc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpicc --with-cxx=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpic++ --with-openmp=0 --download-hypre=yes --download-sowing=yes --download-metis=yes --download-parmetis=yes --download-triangle=yes --download-tetgen=yes --download-ctetgen=yes --download-p4est=yes --download-zlib=yes --download-c2html=yes --download-eigen=yes --download-pragmatic=yes --with-hdf5-dir=/usr/local/Cellar/hdf5/1.10.5_1 --with-cmake-dir=/usr/local/Cellar/cmake/3.15.3[0]PETSC ERROR: #1 PetscObjectGetFortranCallback() line 258 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/sys/objects/inherit.c[0]PETSC ERROR: #2 ourbocofunc() line 141 in /Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC ERROR: #3 DMPlexInsertBoundaryValuesRiemann() line 989 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #4 DMPlexInsertBoundaryValues_Plex() line 1052 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #5 DMPlexInsertBoundaryValues() line 1142 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #6 DMPlexComputeResidual_Internal() line 4524 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #7 DMPlexTSComputeRHSFunctionFVM() line 74 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmplexts.c[0]PETSC ERROR: #8 ourdmtsrhsfunc() line 186 in /Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC ERROR: #9 TSComputeRHSFunction_DMLocal() line 105 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmlocalts.c[0]PETSC ERROR: #10 TSComputeRHSFunction() line 653 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #11 TSSSPStep_RK_3() line 120 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC ERROR: #12 TSStep_SSP() line 208 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC ERROR: #13 TSStep() line 3757 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #14 TSSolve() line 4154 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #15 User provided function() line 0 in User file

Second is about the DMProjectFunction wrapper, that I have done so :

static PetscErrorCode ourdmprojfunc(PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) { PetscObjectUseFortranCallback((DM)ctx, dmprojfunc,

(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*), (&dim,&time,x,&Nf,u,_ctx,&ierr)) } PETSC_EXTERN void dmprojectfunction_(DM *dm, PetscReal *time, void (*func)(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*), void *ctx, InsertMode *mode, Vec X, PetscErrorCode *ierr) { PetscErrorCode (*funcarr[1]) (PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_CLASS, &dmprojfunc, (PetscVoidFunction)func, ctx); funcarr[0] = ourdmprojfunc; *ierr = DMProjectFunction(*dm, *time, funcarr, &ctx, *mode, X); }

This time there is no error because I cannot reach this point in the program, but I am not sure anyways how to write this wrapper, especially because of the double pointers that DMProjectFunction takes as arguments.

Does anyone have any idea what could be going wrong with those two wrappers ?

Thank you very much in advance !!

Thibault

Le ven. 18 déc. 2020 à 11:02, Thibault Bridel-Bertomeu < [email protected]> a écrit :

...

Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry !

The MR will appear in a little while !

Thibault

Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected]> a écrit :

...

Good question. There is a trick to limit the amount of work you need to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected]:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

...

Thibault,

A subdirectory for the example is fine; we have other examples that use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git repository with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec *,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode (*)(void **))PetscViewerAndFormatDestroy); } else { *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx); *ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

With the code in attachement it's better. And I set it up in Fortran as follow :: TS :: timeS TSAdapt :: adapt call TSGetAdapt(timeS, adapt, ierr); CHKERRA(ierr) call MyTSAdaptRegister(timeS, adapt, gamma, cfl, ierr); CHKERRA(ierr) Thibault Le lun. 28 déc. 2020 à 18:02, Thibault Bridel-Bertomeu < [email protected]> a écrit :

Hi Jed,

Thanks for your message. I implemented everything in C as you suggested and it works fine except for one thing : the ts->vec_sol does not seem to get updated when seen from the C code (it is on the Fortran side though the solution is correct). As a result, the time step (that uses among other things the max velocity in the domain) is always at the value it gets from the initial solution. Any idea why ts->vec_sol does not seem to be updated ? (I checked the stepnum and time is updated though , when accessed with TSGetTime for instance).

Cheers, Thibault

Le lun. 28 déc. 2020 à 15:30, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Good morning everyone,

Thank you Barry for the answer, it works now !

I am facing (yet) another situation: the TSAdaptRegister function. In the MR on gitlab, Jed mentioned that sometimes, when function pointers are not stored in PETSc objects, one can use stack memory to pass that pointer from fortran to C.

The issue with stack memory is that when it returns, that memory is invalid. You can't use it in this instance.

I think you're going to have problems implementing a TSAdaptCreate_XYZ in Fortran (because the body of that function will need to access private struct members; see below).

I would implement what you need in C and you can call out to Fortran if you want from inside TSAdaptChoose_YourMethod().

PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt adapt) { TSAdapt_DSP *dsp; PetscErrorCode ierr;

PetscFunctionBegin; ierr = PetscNewLog(adapt,&dsp);CHKERRQ(ierr); adapt->reject_safety = 1.0; /* unused */

adapt->data = (void*)dsp; adapt->ops->choose = TSAdaptChoose_DSP; adapt->ops->setfromoptions = TSAdaptSetFromOptions_DSP; adapt->ops->destroy = TSAdaptDestroy_DSP; adapt->ops->view = TSAdaptView_DSP;

ierr = PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetFilter_C",TSAdaptDSPSetFilter_DSP);CHKERRQ(ierr); ierr = PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetPID_C",TSAdaptDSPSetPID_DSP);CHKERRQ(ierr);

ierr = TSAdaptDSPSetFilter_DSP(adapt,"PI42");CHKERRQ(ierr); ierr = TSAdaptRestart_DSP(adapt);CHKERRQ(ierr); PetscFunctionReturn(0); }

...

Can anyone develop that idea ? Because for TSAdaptRegister, i guess the wrapper would start like :

PETSC_EXTERN void tsadaptregister_(char *sname, void (*func)(TSAdapt*,PetscErrorCode*), PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T snamelen)

but then the C TSAdaptRegister function takes a PetscErrorCode (*func)(TSAdapt) function pointer as argument ... I cannot use any FORTRAN_CALLBACK here since I do not have any object to hook it to, and I could not find a similar situation among the pre-existing wrappers. Does anyone have an idea on how to proceed ?

Thanks !!

Thibault

Le mar. 22 déc. 2020 à 21:20, Barry Smith <[email protected]> a écrit :

...

PetscObjectUseFortranCallback((PetscDS)ctx,

*ierr = PetscObjectSetFortranCallback((PetscObject)*prob

It looks like the problem is that these user provided functions do

not

...

take a PetscDS directly as an argument so the Fortran callback information cannot be obtained from them.

The manual page for PetscDSAddBoundary() says

- ctx - An optional user context for bcFunc

but then when it lists the calling sequence for bcFunc it does not list the ctx as an argument, so either the manual page or code is wrong.

It looks like you make the ctx be the PetscDS prob argument when you call PetscDSAddBoundary

In principle this sounds like it might work. I think you need to track through the debugger to see if the ctx passed to ourbocofunc() is actually the PetscDS prob variable and if not why it is not.

Barry

On Dec 22, 2020, at 5:49 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

I have hit two snags while implementing the missing wrappers necessary to transcribe ex11 to Fortran.

First is about the PetscDSAddBoundary wrapper, that I have done so :

static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, ctx); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, *prob); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); }

but when I call it in the program, with adequate routines, I obtain the following error :

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------[0]PETSC ERROR: Corrupt argument: https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC ERROR: Fortran callback not set on this object[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.[0]PETSC ERROR: Petsc Development GIT revision: v3.14.2-297-gf36a7edeb8 GIT Date: 2020-12-18 04:42:53 +0000[0]PETSC ERROR: ../../../bin/eulerian3D on a named macbook-pro-de-thibault.home by tbridel Sun Dec 20 15:05:15 2020[0]PETSC ERROR: Configure options --with-clean=0 --prefix=/Users/tbridel/Documents/1-CODES/04-PETSC/build --with-make-np=2 --with-windows-graphics=0 --with-debugging=0 --download-fblaslapack --download-mpich-shared=0 --with-x=0 --with-pthread=0 --with-valgrind=0 --PETSC_ARCH=macosx --with-fc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpifort --with-cc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpicc --with-cxx=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpic++ --with-openmp=0 --download-hypre=yes --download-sowing=yes --download-metis=yes --download-parmetis=yes --download-triangle=yes --download-tetgen=yes --download-ctetgen=yes --download-p4est=yes --download-zlib=yes --download-c2html=yes --download-eigen=yes --download-pragmatic=yes --with-hdf5-dir=/usr/local/Cellar/hdf5/1.10.5_1 --with-cmake-dir=/usr/local/Cellar/cmake/3.15.3[0]PETSC ERROR: #1 PetscObjectGetFortranCallback() line 258 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/sys/objects/inherit.c[0]PETSC ERROR: #2 ourbocofunc() line 141 in /Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC ERROR: #3 DMPlexInsertBoundaryValuesRiemann() line 989 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #4 DMPlexInsertBoundaryValues_Plex() line 1052 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #5 DMPlexInsertBoundaryValues() line 1142 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #6 DMPlexComputeResidual_Internal() line 4524 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC ERROR: #7 DMPlexTSComputeRHSFunctionFVM() line 74 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmplexts.c[0]PETSC ERROR: #8 ourdmtsrhsfunc() line 186 in /Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC ERROR: #9 TSComputeRHSFunction_DMLocal() line 105 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmlocalts.c[0]PETSC ERROR: #10 TSComputeRHSFunction() line 653 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #11 TSSSPStep_RK_3() line 120 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC ERROR: #12 TSStep_SSP() line 208 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC ERROR: #13 TSStep() line 3757 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #14 TSSolve() line 4154 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #15 User provided function() line 0 in User file

Second is about the DMProjectFunction wrapper, that I have done so :

static PetscErrorCode ourdmprojfunc(PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) { PetscObjectUseFortranCallback((DM)ctx, dmprojfunc,

(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

(&dim,&time,x,&Nf,u,_ctx,&ierr)) } PETSC_EXTERN void dmprojectfunction_(DM *dm, PetscReal *time, void

(*func)(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

void *ctx, InsertMode *mode, Vec

X,

...

PetscErrorCode *ierr) { PetscErrorCode (*funcarr[1]) (PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_CLASS, &dmprojfunc, (PetscVoidFunction)func, ctx); funcarr[0] = ourdmprojfunc; *ierr = DMProjectFunction(*dm, *time, funcarr, &ctx, *mode, X); }

This time there is no error because I cannot reach this point in the program, but I am not sure anyways how to write this wrapper, especially because of the double pointers that DMProjectFunction takes as arguments.

Does anyone have any idea what could be going wrong with those two wrappers ?

Thank you very much in advance !!

Thibault

Le ven. 18 déc. 2020 à 11:02, Thibault Bridel-Bertomeu < [email protected]> a écrit :

...

Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry !

The MR will appear in a little while !

Thibault

Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected]> a écrit :

...

Good question. There is a trick to limit the amount of work you

need

...

to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected]:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

...

Thibault,

A subdirectory for the example is fine; we have other examples

that

...

use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git

repository

...

with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec

...

*,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void (*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void *mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode

(*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode

...

(*)(void **))PetscViewerAndFormatDestroy); } else { *ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx);

...

*ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx);

...

*ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

I’m using TSSSP steppers. I know the RHS is not zero because when I output the solution from the Fortran code in VTK format, it looks as expected and indeed moves forward with time. I won’t put the code as an attachement there is a little too much but the MR is up to date on the VTK output (at least). Le lun. 28 déc. 2020 à 18:42, Jed Brown <[email protected]> a écrit :

I think I'm not following something. The time stepper implementation (TSStep_XYZ) is where ts->vec_sol should be updated. The controller doesn't have anything to do with it. What TS are you using and how do you know your RHS is nonzero?

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Hi Jed,

Thanks for your message. I implemented everything in C as you suggested and it works fine except for one thing : the ts->vec_sol does not seem to get updated when seen from the C code (it is on the Fortran side though the solution is correct). As a result, the time step (that uses among other things the max velocity in the domain) is always at the value it gets from the initial solution. Any idea why ts->vec_sol does not seem to be updated ? (I checked the stepnum and time is updated though , when accessed with TSGetTime for instance).

Cheers, Thibault

Le lun. 28 déc. 2020 à 15:30, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Good morning everyone,

Thank you Barry for the answer, it works now !

I am facing (yet) another situation: the TSAdaptRegister function. In the MR on gitlab, Jed mentioned that sometimes, when function pointers are not stored in PETSc objects, one can use stack memory to pass that pointer from fortran to C.

The issue with stack memory is that when it returns, that memory is invalid. You can't use it in this instance.

I think you're going to have problems implementing a TSAdaptCreate_XYZ in Fortran (because the body of that function will need to access private struct members; see below).

I would implement what you need in C and you can call out to Fortran if you want from inside TSAdaptChoose_YourMethod().

PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt adapt) { TSAdapt_DSP *dsp; PetscErrorCode ierr;

PetscFunctionBegin; ierr = PetscNewLog(adapt,&dsp);CHKERRQ(ierr); adapt->reject_safety = 1.0; /* unused */

adapt->data = (void*)dsp; adapt->ops->choose = TSAdaptChoose_DSP; adapt->ops->setfromoptions = TSAdaptSetFromOptions_DSP; adapt->ops->destroy = TSAdaptDestroy_DSP; adapt->ops->view = TSAdaptView_DSP;

ierr =

PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetFilter_C",TSAdaptDSPSetFilter_DSP);CHKERRQ(ierr);

...

ierr =

PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetPID_C",TSAdaptDSPSetPID_DSP);CHKERRQ(ierr);

...

ierr = TSAdaptDSPSetFilter_DSP(adapt,"PI42");CHKERRQ(ierr); ierr = TSAdaptRestart_DSP(adapt);CHKERRQ(ierr); PetscFunctionReturn(0); }

...

Can anyone develop that idea ? Because for TSAdaptRegister, i guess

the

...

...

wrapper would start like :

PETSC_EXTERN void tsadaptregister_(char *sname, void (*func)(TSAdapt*,PetscErrorCode*), PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T snamelen)

but then the C TSAdaptRegister function takes a PetscErrorCode (*func)(TSAdapt) function pointer as argument ... I cannot use any FORTRAN_CALLBACK here since I do not have any object to hook it to, and I could not find a similar situation among the pre-existing wrappers. Does anyone have an idea on how to proceed ?

Thanks !!

Thibault

Le mar. 22 déc. 2020 à 21:20, Barry Smith <[email protected]> a écrit :

...

PetscObjectUseFortranCallback((PetscDS)ctx,

*ierr = PetscObjectSetFortranCallback((PetscObject)*prob

It looks like the problem is that these user provided functions do

not

...

take a PetscDS directly as an argument so the Fortran callback information cannot be obtained from them.

The manual page for PetscDSAddBoundary() says

- ctx - An optional user context for bcFunc

but then when it lists the calling sequence for bcFunc it does not list the ctx as an argument, so either the manual page or code is wrong.

It looks like you make the ctx be the PetscDS prob argument when you call PetscDSAddBoundary

In principle this sounds like it might work. I think you need to track through the debugger to see if the ctx passed to ourbocofunc() is actually the PetscDS prob variable and if not why it is not.

Barry

On Dec 22, 2020, at 5:49 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

I have hit two snags while implementing the missing wrappers necessary to transcribe ex11 to Fortran.

First is about the PetscDSAddBoundary wrapper, that I have done so :

static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, ctx); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, *prob); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); }

but when I call it in the program, with adequate routines, I obtain the following error :

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------[0]PETSC ERROR: Corrupt argument: https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC ERROR: Fortran callback not set on this object[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.[0]PETSC ERROR: Petsc Development GIT revision: v3.14.2-297-gf36a7edeb8 GIT Date: 2020-12-18 04:42:53 +0000[0]PETSC ERROR: ../../../bin/eulerian3D on a named macbook-pro-de-thibault.home by tbridel Sun Dec 20 15:05:15 2020[0]PETSC ERROR: Configure options --with-clean=0 --prefix=/Users/tbridel/Documents/1-CODES/04-PETSC/build --with-make-np=2 --with-windows-graphics=0 --with-debugging=0 --download-fblaslapack --download-mpich-shared=0 --with-x=0 --with-pthread=0 --with-valgrind=0 --PETSC_ARCH=macosx --with-fc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpifort --with-cc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpicc --with-cxx=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpic++ --with-openmp=0 --download-hypre=yes --download-sowing=yes --download-metis=yes --download-parmetis=yes --download-triangle=yes --download-tetgen=yes --download-ctetgen=yes --download-p4est=yes --download-zlib=yes --download-c2html=yes --download-eigen=yes --download-pragmatic=yes --with-hdf5-dir=/usr/local/Cellar/hdf5/1.10.5_1 --with-cmake-dir=/usr/local/Cellar/cmake/3.15.3[0]PETSC ERROR: #1 PetscObjectGetFortranCallback() line 258 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/sys/objects/inherit.c[0]PETSC

...

ERROR: #2 ourbocofunc() line 141 in

/Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC

...

ERROR: #3 DMPlexInsertBoundaryValuesRiemann() line 989 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #4 DMPlexInsertBoundaryValues_Plex() line 1052 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #5 DMPlexInsertBoundaryValues() line 1142 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #6 DMPlexComputeResidual_Internal() line 4524 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #7 DMPlexTSComputeRHSFunctionFVM() line 74 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmplexts.c[0]PETSC

...

ERROR: #8 ourdmtsrhsfunc() line 186 in

/Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC

...

ERROR: #9 TSComputeRHSFunction_DMLocal() line 105 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmlocalts.c[0]PETSC

...

ERROR: #10 TSComputeRHSFunction() line 653 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #11 TSSSPStep_RK_3() line 120 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC

...

ERROR: #12 TSStep_SSP() line 208 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC

...

ERROR: #13 TSStep() line 3757 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #14 TSSolve() line 4154 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #15 User provided function() line 0 in User file

...

...

Second is about the DMProjectFunction wrapper, that I have done so :

static PetscErrorCode ourdmprojfunc(PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) { PetscObjectUseFortranCallback((DM)ctx, dmprojfunc,

(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

...

...

(&dim,&time,x,&Nf,u,_ctx,&ierr)) } PETSC_EXTERN void dmprojectfunction_(DM *dm, PetscReal *time, void

(*func)(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

...

...

void *ctx, InsertMode *mode,

Vec X,

...

PetscErrorCode *ierr) { PetscErrorCode (*funcarr[1]) (PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_CLASS, &dmprojfunc, (PetscVoidFunction)func, ctx); funcarr[0] = ourdmprojfunc; *ierr = DMProjectFunction(*dm, *time, funcarr, &ctx, *mode, X); }

This time there is no error because I cannot reach this point in the program, but I am not sure anyways how to write this wrapper, especially because of the double pointers that DMProjectFunction takes as arguments.

Does anyone have any idea what could be going wrong with those two wrappers ?

Thank you very much in advance !!

Thibault

Le ven. 18 déc. 2020 à 11:02, Thibault Bridel-Bertomeu < [email protected]> a écrit :

...

Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry !

The MR will appear in a little while !

Thibault

Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected]> a écrit :

...

Good question. There is a trick to limit the amount of work you

need

...

to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected]:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

...

Thibault,

A subdirectory for the example is fine; we have other examples

that

...

use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git

repository

...

with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec

...

...

...

...

...

...

...

...

*,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void

(*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void

...

*mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode

(*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode

...

...

...

...

...

...

...

...

(*)(void **))PetscViewerAndFormatDestroy); } else { *ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx);

...

...

...

...

...

...

...

...

*ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx);

...

...

...

...

...

...

...

...

*ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

I'd need some form of reproducible test case to debug. It sounds like you're viewing different vectors, but I probably won't be able to track it down without a reproducible test case. A gdb session might be most helpful to identify which vectors are where when you get the unexpected condition. Thibault Bridel-Bertomeu <[email protected]> writes:

However in Fortran, for the VTK output, it is true I do not view ts->vec_sol directly but rather a Vec called sol that I use in the TSSolve(ts, sol) call ... but I suppose it should be equivalent shouldn’t it ?

Le lun. 28 déc. 2020 à 18:42, Jed Brown <[email protected]> a écrit :

...

I think I'm not following something. The time stepper implementation (TSStep_XYZ) is where ts->vec_sol should be updated. The controller doesn't have anything to do with it. What TS are you using and how do you know your RHS is nonzero?

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Hi Jed,

Thanks for your message. I implemented everything in C as you suggested and it works fine except for one thing : the ts->vec_sol does not seem to get updated when seen from the C code (it is on the Fortran side though the solution is correct). As a result, the time step (that uses among other things the max velocity in the domain) is always at the value it gets from the initial solution. Any idea why ts->vec_sol does not seem to be updated ? (I checked the stepnum and time is updated though , when accessed with TSGetTime for instance).

Cheers, Thibault

Le lun. 28 déc. 2020 à 15:30, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Good morning everyone,

Thank you Barry for the answer, it works now !

I am facing (yet) another situation: the TSAdaptRegister function. In the MR on gitlab, Jed mentioned that sometimes, when function pointers are not stored in PETSc objects, one can use stack memory to pass that pointer from fortran to C.

The issue with stack memory is that when it returns, that memory is invalid. You can't use it in this instance.

I think you're going to have problems implementing a TSAdaptCreate_XYZ in Fortran (because the body of that function will need to access private struct members; see below).

I would implement what you need in C and you can call out to Fortran if you want from inside TSAdaptChoose_YourMethod().

PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt adapt) { TSAdapt_DSP *dsp; PetscErrorCode ierr;

PetscFunctionBegin; ierr = PetscNewLog(adapt,&dsp);CHKERRQ(ierr); adapt->reject_safety = 1.0; /* unused */

adapt->data = (void*)dsp; adapt->ops->choose = TSAdaptChoose_DSP; adapt->ops->setfromoptions = TSAdaptSetFromOptions_DSP; adapt->ops->destroy = TSAdaptDestroy_DSP; adapt->ops->view = TSAdaptView_DSP;

ierr =

PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetFilter_C",TSAdaptDSPSetFilter_DSP);CHKERRQ(ierr);

...

ierr =

PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetPID_C",TSAdaptDSPSetPID_DSP);CHKERRQ(ierr);

...

ierr = TSAdaptDSPSetFilter_DSP(adapt,"PI42");CHKERRQ(ierr); ierr = TSAdaptRestart_DSP(adapt);CHKERRQ(ierr); PetscFunctionReturn(0); }

...

Can anyone develop that idea ? Because for TSAdaptRegister, i guess

the

...

...

wrapper would start like :

PETSC_EXTERN void tsadaptregister_(char *sname, void (*func)(TSAdapt*,PetscErrorCode*), PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T snamelen)

but then the C TSAdaptRegister function takes a PetscErrorCode (*func)(TSAdapt) function pointer as argument ... I cannot use any FORTRAN_CALLBACK here since I do not have any object to hook it to, and I could not find a similar situation among the pre-existing wrappers. Does anyone have an idea on how to proceed ?

Thanks !!

Thibault

Le mar. 22 déc. 2020 à 21:20, Barry Smith <[email protected]> a écrit :

...

PetscObjectUseFortranCallback((PetscDS)ctx,

*ierr = PetscObjectSetFortranCallback((PetscObject)*prob

It looks like the problem is that these user provided functions do

not

...

take a PetscDS directly as an argument so the Fortran callback information cannot be obtained from them.

The manual page for PetscDSAddBoundary() says

- ctx - An optional user context for bcFunc

but then when it lists the calling sequence for bcFunc it does not list the ctx as an argument, so either the manual page or code is wrong.

It looks like you make the ctx be the PetscDS prob argument when you call PetscDSAddBoundary

In principle this sounds like it might work. I think you need to track through the debugger to see if the ctx passed to ourbocofunc() is actually the PetscDS prob variable and if not why it is not.

Barry

On Dec 22, 2020, at 5:49 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

I have hit two snags while implementing the missing wrappers necessary to transcribe ex11 to Fortran.

First is about the PetscDSAddBoundary wrapper, that I have done so :

static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, ctx); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, *prob); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); }

but when I call it in the program, with adequate routines, I obtain the following error :

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------[0]PETSC ERROR: Corrupt argument: https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC ERROR: Fortran callback not set on this object[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.[0]PETSC ERROR: Petsc Development GIT revision: v3.14.2-297-gf36a7edeb8 GIT Date: 2020-12-18 04:42:53 +0000[0]PETSC ERROR: ../../../bin/eulerian3D on a named macbook-pro-de-thibault.home by tbridel Sun Dec 20 15:05:15 2020[0]PETSC ERROR: Configure options --with-clean=0 --prefix=/Users/tbridel/Documents/1-CODES/04-PETSC/build --with-make-np=2 --with-windows-graphics=0 --with-debugging=0 --download-fblaslapack --download-mpich-shared=0 --with-x=0 --with-pthread=0 --with-valgrind=0 --PETSC_ARCH=macosx --with-fc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpifort --with-cc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpicc --with-cxx=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpic++ --with-openmp=0 --download-hypre=yes --download-sowing=yes --download-metis=yes --download-parmetis=yes --download-triangle=yes --download-tetgen=yes --download-ctetgen=yes --download-p4est=yes --download-zlib=yes --download-c2html=yes --download-eigen=yes --download-pragmatic=yes --with-hdf5-dir=/usr/local/Cellar/hdf5/1.10.5_1 --with-cmake-dir=/usr/local/Cellar/cmake/3.15.3[0]PETSC ERROR: #1 PetscObjectGetFortranCallback() line 258 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/sys/objects/inherit.c[0]PETSC

...

ERROR: #2 ourbocofunc() line 141 in

/Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC

...

ERROR: #3 DMPlexInsertBoundaryValuesRiemann() line 989 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #4 DMPlexInsertBoundaryValues_Plex() line 1052 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #5 DMPlexInsertBoundaryValues() line 1142 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #6 DMPlexComputeResidual_Internal() line 4524 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #7 DMPlexTSComputeRHSFunctionFVM() line 74 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmplexts.c[0]PETSC

...

ERROR: #8 ourdmtsrhsfunc() line 186 in

/Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC

...

ERROR: #9 TSComputeRHSFunction_DMLocal() line 105 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmlocalts.c[0]PETSC

...

ERROR: #10 TSComputeRHSFunction() line 653 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #11 TSSSPStep_RK_3() line 120 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC

...

ERROR: #12 TSStep_SSP() line 208 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC

...

ERROR: #13 TSStep() line 3757 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #14 TSSolve() line 4154 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #15 User provided function() line 0 in User file

...

...

Second is about the DMProjectFunction wrapper, that I have done so :

static PetscErrorCode ourdmprojfunc(PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) { PetscObjectUseFortranCallback((DM)ctx, dmprojfunc,

(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

...

...

(&dim,&time,x,&Nf,u,_ctx,&ierr)) } PETSC_EXTERN void dmprojectfunction_(DM *dm, PetscReal *time, void

(*func)(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

...

...

void *ctx, InsertMode *mode,

Vec X,

...

PetscErrorCode *ierr) { PetscErrorCode (*funcarr[1]) (PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_CLASS, &dmprojfunc, (PetscVoidFunction)func, ctx); funcarr[0] = ourdmprojfunc; *ierr = DMProjectFunction(*dm, *time, funcarr, &ctx, *mode, X); }

This time there is no error because I cannot reach this point in the program, but I am not sure anyways how to write this wrapper, especially because of the double pointers that DMProjectFunction takes as arguments.

Does anyone have any idea what could be going wrong with those two wrappers ?

Thank you very much in advance !!

Thibault

Le ven. 18 déc. 2020 à 11:02, Thibault Bridel-Bertomeu < [email protected]> a écrit :

...

Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry !

The MR will appear in a little while !

Thibault

Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected]> a écrit :

...

Good question. There is a trick to limit the amount of work you

need

...

to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected]:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

...

Thibault,

A subdirectory for the example is fine; we have other examples

that

...

use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git

repository

...

with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec

...

...

...

...

...

...

...

...

*,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void

(*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void

...

*mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode

(*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode

...

...

...

...

...

...

...

...

(*)(void **))PetscViewerAndFormatDestroy); } else { *ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx);

...

...

...

...

...

...

...

...

*ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx);

...

...

...

...

...

...

...

...

*ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

Oh okay, from a monitor function it always uses the latest up to date solution ! Then my mistake, when I visualize the VTK it comes from a monitor so it must be the latest solution and not the Vec sol as I was saying, sorry. I’ll check out why in the C part TSGetSolutuon always yields the initial vector rather than the current one ! Keep you in the loop ! Thankss !! Le mar. 29 déc. 2020 à 04:39, Barry Smith <[email protected]> a écrit :

On Dec 28, 2020, at 11:52 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

However in Fortran, for the VTK output, it is true I do not view ts->vec_sol directly but rather a Vec called sol that I use in the TSSolve(ts, sol) call ... but I suppose it should be equivalent shouldn’t it ?

This can be tricky. As the integrator runs it may not always be updating the vector you pass into TSSolve. It is only when TSSolve returns that the values in sol are guaranteed to match the TS solution. It is better if your VTK output code call something to get get the current solution and not use your Fortran sol vector.

For example if you register a monitor function one of its arguments is always the correct current solution. I think you can also call TSGetSolution to return a vector that contains the current solution.

Barry

Le lun. 28 déc. 2020 à 18:42, Jed Brown <[email protected]> a écrit :

...

I think I'm not following something. The time stepper implementation (TSStep_XYZ) is where ts->vec_sol should be updated. The controller doesn't have anything to do with it. What TS are you using and how do you know your RHS is nonzero?

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Hi Jed,

Thanks for your message. I implemented everything in C as you suggested and it works fine except for one thing : the ts->vec_sol does not seem to get updated when seen from the C code (it is on the Fortran side though the solution is correct). As a result, the time step (that uses among other things the max velocity in the domain) is always at the value it gets from the initial solution. Any idea why ts->vec_sol does not seem to be updated ? (I checked the stepnum and time is updated though , when accessed with TSGetTime for instance).

Cheers, Thibault

Le lun. 28 déc. 2020 à 15:30, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Good morning everyone,

Thank you Barry for the answer, it works now !

I am facing (yet) another situation: the TSAdaptRegister function. In the MR on gitlab, Jed mentioned that sometimes, when function pointers are not stored in PETSc objects, one can use stack memory to pass that pointer from fortran to C.

The issue with stack memory is that when it returns, that memory is invalid. You can't use it in this instance.

I think you're going to have problems implementing a TSAdaptCreate_XYZ in Fortran (because the body of that function will need to access private struct members; see below).

I would implement what you need in C and you can call out to Fortran if you want from inside TSAdaptChoose_YourMethod().

PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt adapt) { TSAdapt_DSP *dsp; PetscErrorCode ierr;

PetscFunctionBegin; ierr = PetscNewLog(adapt,&dsp);CHKERRQ(ierr); adapt->reject_safety = 1.0; /* unused */

adapt->data = (void*)dsp; adapt->ops->choose = TSAdaptChoose_DSP; adapt->ops->setfromoptions = TSAdaptSetFromOptions_DSP; adapt->ops->destroy = TSAdaptDestroy_DSP; adapt->ops->view = TSAdaptView_DSP;

ierr =

PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetFilter_C",TSAdaptDSPSetFilter_DSP);CHKERRQ(ierr);

...

ierr =

PetscObjectComposeFunction((PetscObject)adapt,"TSAdaptDSPSetPID_C",TSAdaptDSPSetPID_DSP);CHKERRQ(ierr);

...

ierr = TSAdaptDSPSetFilter_DSP(adapt,"PI42");CHKERRQ(ierr); ierr = TSAdaptRestart_DSP(adapt);CHKERRQ(ierr); PetscFunctionReturn(0); }

...

Can anyone develop that idea ? Because for TSAdaptRegister, i guess

the

...

...

wrapper would start like :

PETSC_EXTERN void tsadaptregister_(char *sname, void (*func)(TSAdapt*,PetscErrorCode*), PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T snamelen)

but then the C TSAdaptRegister function takes a PetscErrorCode (*func)(TSAdapt) function pointer as argument ... I cannot use any FORTRAN_CALLBACK here since I do not have any object to hook it to, and I could not find a similar situation among the pre-existing wrappers. Does anyone have an idea on how to proceed ?

Thanks !!

Thibault

Le mar. 22 déc. 2020 à 21:20, Barry Smith <[email protected]> a écrit :

...

PetscObjectUseFortranCallback((PetscDS)ctx,

*ierr = PetscObjectSetFortranCallback((PetscObject)*prob

It looks like the problem is that these user provided functions

do not

...

take a PetscDS directly as an argument so the Fortran callback information cannot be obtained from them.

The manual page for PetscDSAddBoundary() says

- ctx - An optional user context for bcFunc

but then when it lists the calling sequence for bcFunc it does not list the ctx as an argument, so either the manual page or code is wrong.

It looks like you make the ctx be the PetscDS prob argument when you call PetscDSAddBoundary

In principle this sounds like it might work. I think you need to track through the debugger to see if the ctx passed to ourbocofunc() is actually the PetscDS prob variable and if not why it is not.

Barry

On Dec 22, 2020, at 5:49 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

I have hit two snags while implementing the missing wrappers necessary to transcribe ex11 to Fortran.

First is about the PetscDSAddBoundary wrapper, that I have done so :

static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*,const PetscReal*,const PetscReal*,const PetscScalar*,const PetscScalar*,void*,PetscErrorCode*), (&time,c,n,a_xI,a_xG,ctx,&ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, ctx); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, *prob); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); }

but when I call it in the program, with adequate routines, I obtain the following error :

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------[0]PETSC ERROR: Corrupt argument: https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC ERROR: Fortran callback not set on this object[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.[0]PETSC ERROR: Petsc Development GIT revision: v3.14.2-297-gf36a7edeb8 GIT Date: 2020-12-18 04:42:53 +0000[0]PETSC ERROR: ../../../bin/eulerian3D on a named macbook-pro-de-thibault.home by tbridel Sun Dec 20 15:05:15 2020[0]PETSC ERROR: Configure options --with-clean=0 --prefix=/Users/tbridel/Documents/1-CODES/04-PETSC/build --with-make-np=2 --with-windows-graphics=0 --with-debugging=0 --download-fblaslapack --download-mpich-shared=0 --with-x=0 --with-pthread=0 --with-valgrind=0 --PETSC_ARCH=macosx --with-fc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpifort --with-cc=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpicc --with-cxx=/usr/local/Cellar/open-mpi/4.0.1_2/bin/mpic++ --with-openmp=0 --download-hypre=yes --download-sowing=yes --download-metis=yes --download-parmetis=yes --download-triangle=yes --download-tetgen=yes --download-ctetgen=yes --download-p4est=yes --download-zlib=yes --download-c2html=yes --download-eigen=yes --download-pragmatic=yes --with-hdf5-dir=/usr/local/Cellar/hdf5/1.10.5_1 --with-cmake-dir=/usr/local/Cellar/cmake/3.15.3[0]PETSC ERROR: #1 PetscObjectGetFortranCallback() line 258 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/sys/objects/inherit.c[0]PETSC

...

ERROR: #2 ourbocofunc() line 141 in

/Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC

...

ERROR: #3 DMPlexInsertBoundaryValuesRiemann() line 989 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #4 DMPlexInsertBoundaryValues_Plex() line 1052 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #5 DMPlexInsertBoundaryValues() line 1142 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #6 DMPlexComputeResidual_Internal() line 4524 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/dm/impls/plex/plexfem.c[0]PETSC

...

ERROR: #7 DMPlexTSComputeRHSFunctionFVM() line 74 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmplexts.c[0]PETSC

...

ERROR: #8 ourdmtsrhsfunc() line 186 in

/Users/tbridel/Documents/1-CODES/59-EULERIAN3D/sources/petsc_wrapping/wrapper_petsc.c[0]PETSC

...

ERROR: #9 TSComputeRHSFunction_DMLocal() line 105 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/utils/dmlocalts.c[0]PETSC

...

ERROR: #10 TSComputeRHSFunction() line 653 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #11 TSSSPStep_RK_3() line 120 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC

...

ERROR: #12 TSStep_SSP() line 208 in

/Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/impls/explicit/ssp/ssp.c[0]PETSC

...

ERROR: #13 TSStep() line 3757 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #14 TSSolve() line 4154 in /Users/tbridel/Documents/1-CODES/04-PETSC/src/ts/interface/ts.c[0]PETSC ERROR: #15 User provided function() line 0 in User file

...

...

Second is about the DMProjectFunction wrapper, that I have done so :

static PetscErrorCode ourdmprojfunc(PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) { PetscObjectUseFortranCallback((DM)ctx, dmprojfunc,

(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

...

...

(&dim,&time,x,&Nf,u,_ctx,&ierr)) } PETSC_EXTERN void dmprojectfunction_(DM *dm, PetscReal *time, void

(*func)(PetscInt*,PetscReal*,PetscReal*,PetscInt*,PetscScalar*,void*,PetscErrorCode*),

...

...

...

void *ctx, InsertMode *mode,

Vec X,

...

PetscErrorCode *ierr) { PetscErrorCode (*funcarr[1]) (PetscInt dim, PetscReal time, PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); *ierr = PetscObjectSetFortranCallback((PetscObject)*dm, PETSC_FORTRAN_CALLBACK_CLASS, &dmprojfunc, (PetscVoidFunction)func, ctx); funcarr[0] = ourdmprojfunc; *ierr = DMProjectFunction(*dm, *time, funcarr, &ctx, *mode, X); }

This time there is no error because I cannot reach this point in the program, but I am not sure anyways how to write this wrapper, especially because of the double pointers that DMProjectFunction takes as arguments.

Does anyone have any idea what could be going wrong with those two wrappers ?

Thank you very much in advance !!

Thibault

Le ven. 18 déc. 2020 à 11:02, Thibault Bridel-Bertomeu < [email protected]> a écrit :

...

Aah that is a nice trick, I was getting ready to fork, clone the fork and redo the work, but that worked fine ! Thank you Barry !

The MR will appear in a little while !

Thibault

Le ven. 18 déc. 2020 à 10:16, Barry Smith <[email protected]> a écrit :

...

Good question. There is a trick to limit the amount of work you

need

...

to do with a new fork after you have already made changes with a PETSc clone, but it looks like we do not document this clearly in the webpages. (I couldn't find it).

Yes, you do need to make a fork, but after you have made the fork on the GitLab website (and have done nothing on your machine) edit the file $PETSC_DIR/.git/config for your clone on your machine

Locate the line that has url = [email protected]:petsc/petsc.git (this may have an https at the beginning of the line)

Change this line to point to the fork url instead with git@ not https, which will be pretty much the same URL but with your user id instead of petsc in the address. Then git push and it will push to your fork.

Now you changes will be in your fork and you can make the MR from your fork URL on Gitlab. (In other words this editing trick converts your PETSc clone on your machine to a PETSc fork).

I hope I have explained this clearly enough it goes smoothly.

Barry

On Dec 18, 2020, at 3:00 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello Barry,

I'll start the MR as soon as possible then so that specialists can indeed have a look. Do I have to fork PETSc to start a MR or are PETSc repo settings such that can I push a branch from the PETSc clone I got ?

Thibault

Le mer. 16 déc. 2020 à 07:47, Barry Smith <[email protected]> a écrit :

...

Thibault,

A subdirectory for the example is fine; we have other examples

that

...

use subdirectories and multiple files.

Note: even if you don't have something completely working you can still make MR and list it as DRAFT request for comments, some other PETSc members who understand the packages you are using and Fortran better than I may be able to help as you develop the code.

Barry

On Dec 16, 2020, at 12:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

Thank you Barry for the feedback. OK, yes I'll work up an MR as soon as I have got something working. By the way, does the fortran-version of the example have to be a single file ? If my push contains a directory with several files (different modules and the main), and the Makefile that goes with it, is that ok ?

Thibault Bridel-Bertomeu

Le mer. 16 déc. 2020 à 04:46, Barry Smith <[email protected]> a écrit :

...

This is great. If you make a branch off of the PETSc git

repository

...

with these additions and work on ex11 you can make a merge request and we can run the code easily on all our test systems (for security reasons one of use needs to launch the tests from your MR). https://docs.petsc.org/en/latest/developers/integration/

Barry

On Dec 15, 2020, at 5:35 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Hello everyone,

So far, I have the wrappers in the files attached to this e-mail. I still do not know if they work properly - at least the code compiles and the calls to the wrapped-subroutine do not fail - but I wanted to put this here in case someone sees something really wrong with it already.

Thank you again for your help, I'll try to post updates of the F90 version of ex11 regularly in this thread.

Stay safe,

Thibault Bridel-Bertomeu

Le dim. 13 déc. 2020 à 16:39, Jed Brown <[email protected]> a écrit :

...

Thibault Bridel-Bertomeu <[email protected]> writes:

...

Thank you Mark for your answer.

I am not sure what you think could be in the setBC1 routine ? How to make the connection with the PetscDS ?

On the other hand, I actually found after a while TSMonitorSet has a fortran wrapper, and it does take as arguments two function pointers, so I guess it is possible ? Although I am not sure exactly how to play with the PetscObjectSetFortranCallback & PetscObjectUseFortranCallback macros - could anybody advise please ?

tsmonitorset_ is a good example to follow. In your file, create one of these static structs with a member for each callback. These are IDs that will be used as keys for Fortran callbacks and their contexts. The salient parts of the file are below.

static struct { PetscFortranCallbackId prestep; PetscFortranCallbackId poststep; PetscFortranCallbackId rhsfunction; PetscFortranCallbackId rhsjacobian; PetscFortranCallbackId ifunction; PetscFortranCallbackId ijacobian; PetscFortranCallbackId monitor; PetscFortranCallbackId mondestroy; PetscFortranCallbackId transform; #if defined(PETSC_HAVE_F90_2PTR_ARG) PetscFortranCallbackId function_pgiptr; #endif } _cb;

/* Note ctx is the same as ts so we need to get the Fortran context out of the TS; this gets put in _ctx using the callback ID */ static PetscErrorCode ourmonitor(TS ts,PetscInt i,PetscReal d,Vec v,void *ctx) {

PetscObjectUseFortranCallback(ts,_cb.monitor,(TS*,PetscInt*,PetscReal*,Vec

...

...

...

...

...

...

...

...

*,void*,PetscErrorCode*),(&ts,&i,&d,&v,_ctx,&ierr)); }

Then follow as in tsmonitorset_, which sets two callbacks.

PETSC_EXTERN void tsmonitorset_(TS *ts,void

(*func)(TS*,PetscInt*,PetscReal*,Vec*,void*,PetscErrorCode*),void

...

*mctx,void (*d)(void*,PetscErrorCode*),PetscErrorCode *ierr) { CHKFORTRANNULLFUNCTION(d); if ((PetscVoidFunction)func == (PetscVoidFunction) tsmonitordefault_) { *ierr = TSMonitorSet(*ts,(PetscErrorCode

(*)(TS,PetscInt,PetscReal,Vec,void*))TSMonitorDefault,*(PetscViewerAndFormat**)mctx,(PetscErrorCode

...

...

...

...

...

...

...

...

(*)(void **))PetscViewerAndFormatDestroy); } else { *ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.monitor,(PetscVoidFunction)func,mctx);

...

...

...

...

...

...

...

...

*ierr =

PetscObjectSetFortranCallback((PetscObject)*ts,PETSC_FORTRAN_CALLBACK_CLASS,&_cb.mondestroy,(PetscVoidFunction)d,mctx);

...

...

...

...

...

...

...

...

*ierr = TSMonitorSet(*ts,ourmonitor,*ts,ourmonitordestroy); } }

<wrapper_petsc.h90><wrapper_petsc.c>

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

-- Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

--

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

Hello Matthew, Thank you for your answer and the pointer to SNESSetFunction(). Is the PETSC_F90_2PTR_PROTO(ptr) necessary for every wrapper of that category ? So far, for the PetscDSSetRiemannSolver, I did this : #ifdef PETSC_HAVE_FORTRAN_CAPS #define petscdssetriemannsolver_ PETSCDSSETRIEMANNSOLVER #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) && !defined(FORTRANDOUBLEUNDERSCORE) #define petscdssetriemannsolver_ petscdssetriemannsolver #endif PetscFortranCallbackId riemannsolver; static PetscErrorCode ourriemannsolver(PetscInt dim, PetscInt Nf, PetscReal x[], PetscReal n[], PetscScalar uL[], PetscScalar uR[], PetscInt numConstants, PetscScalar constants[], PetscScalar flux[], void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, riemannsolver, (PetscInt*, PetscInt*, PetscReal*, PetscReal*, PetscScalar*, PetscScalar*, PetscInt*, PetscScalar*, PetscScalar*, void*, PetscErrorCode*), (&dim, &Nf, x, n, uL, uR, &numConstants, constants, flux, _ctx, &ierr)); } PETSC_EXTERN void petscdssetriemannsolver_(PetscDS *prob, PetscInt *f, void (*rs)(PetscInt *dim, PetscInt *Nf, PetscReal x[], PetscReal n[], PetscScalar uL[], PetscScalar uR[], PetscInt *numConstants, PetscScalar constants[], PetscScalar flux[], void *ctx, PetscErrorCode *jerr), PetscErrorCode *ierr) { *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &riemannsolver, (PetscVoidFunction)rs, prob); *ierr = PetscDSSetRiemannSolver(*prob, *f, (void*)ourriemannsolver); } It compiles, but I have no idea whether it works or not, and I won't have until the program is entirely built (with the TS and all) and it actually has to call a Riemann solver added with that routine. What do you think ? Thibault Le dim. 13 déc. 2020 à 15:29, Matthew Knepley <[email protected]> a écrit :

On Sun, Dec 13, 2020 at 9:17 AM Mark Adams <[email protected]> wrote:

...

I don't think function pointers to PETSc (DS and DM) methods are not going to work in Fortran:

ierr = PetscDSAddBoundary(prob, DM_BC_NATURAL_RIEMANN, "inflow", "Face Sets", 0, 0, NULL, (void (*)(void)) PhysicsBoundary_Advect_Inflow, NULL, ALEN(inflowids), inflowids, phys);CHKERRQ(ierr);

You could write a funcs.c file that you call from your fortran code, like, call setBC1(prob,...,ierr)

and put PhysicsBoundary_Advect_Inflow and setBC1 in funcs.c, for instance.

Wrappers for functions that that function arguments are possible, but more involved. We have to create an internal struct that holds Fortran function pointers, and then call them with the right arguments. This has been done for SNESSetFunction() for instance, so we would start emulating that.

Thanks,

Matt

...

On Sun, Dec 13, 2020 at 5:32 AM Thibault Bridel-Bertomeu < [email protected]> wrote:

...

Good morning all,

Thank you Barry for your answer. I started adding some interfaces (PetscFVSetComponentName, PetscFVView & PetscFVSetType so far) and I have to say I think it is working quite well, but the prototypes of those functions are still quite "simple". I am stuck at how to implement the wrappers for PetscDSSetRiemannSolver and PetscDSSetContext though, especially on how to pass a function as an argument to PetscDSSetRiemannSolver ... Are there any similar functions that may already have their wrappers ?

Thank you very much,

Thibault

Le sam. 12 déc. 2020 à 23:28, Barry Smith <[email protected]> a écrit :

...

On Dec 12, 2020, at 2:59 PM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear Jed, dear Barry,

Thank you for the fast answers !

If I have any success I will make sure to make a pull request to provide some version of ex11 in Fortran.

Regarding the stubs, I admit I started looking in that direction to add the missing wrappers but I am not sure I fully understand the process yet. For each C function, I gotta provide a Fortran interface in a .h90 file as well as a C function that has a Fortran-like prototype and calls the C function - right ?

Yes,

However there are a few things I could not find / understand yet. For instance, it appears that for C functions that have character string arguments take an extra argument in their Fortran-like-prototype-wrapper, namely the length of the string. Is that passed automatically ? I couldn’t find where it could come from ...

This secret argument is put automatically by the Fortran compiler.

Another thing is for functions like PetscFVView. I guess the wrapping is less straightforward because I tried a quick something and it segfault’ed. I couldnt find the wrapper for DMView although there is such a routine in Fortran too. Could you please detail how to wrap such functions ?

PETSC_EXTERN void dmview_(DM *da,PetscViewer *vin,PetscErrorCode *ierr) { PetscViewer v; PetscPatchDefaultViewers_Fortran(vin,v); *ierr = DMView(*da,v); }

dm/interface/ftn-custom/zdmf.c

Thank you very much again,

Thibault Bridel-Bertomeu

Le sam. 12 déc. 2020 à 21:48, Barry Smith <[email protected]> a écrit :

...

PETSc Fortran interfaces are a combination of automatically generated and manually generated.

For any C PETSc function if the manual page begins with /*@ it generates the Fortran interface automatically (make allfortranstubs). If it begins /*@C then either the Fortran interface is done manually or is missing.

C functions that have character string arguments or function arguments (or a few other special cases) need to be manually provided. The automatically generated stubs go in the directory ftn-auto while manually generated ones go in the directory fin-custom.

Perhaps you could first generate a list of "missing" Fortran stubs and then for each stub determine why it is missing and if it can be provided. Some are likely easy to provide but a few (involving function arguments) will be more involved. Once you have all the stubs available translating ex11.c becomes straightforward.

Barry

On Dec 12, 2020, at 9:30 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

Is there somewhere a version of the TS tutorial ex11.c in Fortran ? I am looking into building in F90 (let's say that it is an unavoidable constraint) an unstructured 3D solver of the Euler equations using the "new" features of PETSc - mostly DMPlex & PetscFV - but I think there are some interfaces missing and I find it hard to find workarounds in Fortran. I would be grateful if anyone could please give me some pointers ...

Thank you very much in advance,

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

--

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>

OK thank you ! For the PetscDSAddBoundary I am not sure either ... I did this for now : #ifdef PETSC_HAVE_FORTRAN_CAPS #define petscdsaddboundary_ PETSCDSADDBOUNDARY #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) && !defined(FORTRANDOUBLEUNDERSCORE) #define petscdsaddboundary_ petscdsaddboundary #endif PetscFortranCallbackId bocofunc, bocofunc_time; static PetscErrorCode ourbocofunc(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc, (PetscReal*, const PetscReal*, const PetscReal*, const PetscScalar*, const PetscScalar*, void*, PetscErrorCode*), (&time, c, n, a_xI, a_xG, ctx, &ierr)); } static PetscErrorCode ourbocofunc_time(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *a_xI, const PetscScalar *a_xG, void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, bocofunc_time, (PetscReal*, const PetscReal*, const PetscReal*, const PetscScalar*, const PetscScalar*, void*, PetscErrorCode*), (&time, c, n, a_xI, a_xG, ctx, &ierr)); } PETSC_EXTERN void petscdsaddboundary_(PetscDS *prob, DMBoundaryConditionType *type, char *name, char *labelname, PetscInt *field, PetscInt *numcomps, PetscInt *comps, void (*bcFunc)(void), void (*bcFunc_t)(void), PetscInt *numids, const PetscInt *ids, void *ctx, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T namelen, PETSC_FORTRAN_CHARLEN_T labelnamelen) { char *newname, *newlabelname; FIXCHAR(name, namelen, newname); FIXCHAR(labelname, labelnamelen, newlabelname); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc, (PetscVoidFunction)bcFunc, prob); *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &bocofunc_time, (PetscVoidFunction)bcFunc_t, prob); *ierr = PetscDSAddBoundary(*prob, *type, newname, newlabelname, *field, *numcomps, comps, (void (*)(void))ourbocofunc, (void (*)(void))ourbocofunc_time, *numids, ids, ctx); FREECHAR(name, newname); FREECHAR(labelname, newlabelname); } But I do not know how to handle the two char* is the argument list : are there going to be two PETSC_FORTRAN_CHARLEN_T like I wrote ? Thanks again, Thibault Le dim. 13 déc. 2020 à 15:52, Matthew Knepley <[email protected]> a écrit :

On Sun, Dec 13, 2020 at 9:35 AM Thibault Bridel-Bertomeu < [email protected]> wrote:

...

Hello Matthew,

Thank you for your answer and the pointer to SNESSetFunction(). Is the PETSC_F90_2PTR_PROTO(ptr) necessary for every wrapper of that category ?

So far, for the PetscDSSetRiemannSolver, I did this :

#ifdef PETSC_HAVE_FORTRAN_CAPS #define petscdssetriemannsolver_ PETSCDSSETRIEMANNSOLVER #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) && !defined(FORTRANDOUBLEUNDERSCORE) #define petscdssetriemannsolver_ petscdssetriemannsolver #endif

PetscFortranCallbackId riemannsolver;

static PetscErrorCode ourriemannsolver(PetscInt dim, PetscInt Nf, PetscReal x[], PetscReal n[], PetscScalar uL[], PetscScalar uR[], PetscInt numConstants, PetscScalar constants[], PetscScalar flux[], void *ctx) { PetscObjectUseFortranCallback((PetscDS)ctx, riemannsolver, (PetscInt*, PetscInt*, PetscReal*, PetscReal*, PetscScalar*, PetscScalar*, PetscInt*, PetscScalar*, PetscScalar*, void*, PetscErrorCode*), (&dim, &Nf, x, n, uL, uR, &numConstants, constants, flux, _ctx, &ierr)); }

PETSC_EXTERN void petscdssetriemannsolver_(PetscDS *prob, PetscInt *f, void (*rs)(PetscInt *dim, PetscInt *Nf, PetscReal x[], PetscReal n[], PetscScalar uL[], PetscScalar uR[], PetscInt *numConstants, PetscScalar constants[], PetscScalar flux[], void *ctx, PetscErrorCode *jerr), PetscErrorCode *ierr) { *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &riemannsolver, (PetscVoidFunction)rs, prob); *ierr = PetscDSSetRiemannSolver(*prob, *f, (void*)ourriemannsolver); }

It compiles, but I have no idea whether it works or not, and I won't have until the program is entirely built (with the TS and all) and it actually has to call a Riemann solver added with that routine. What do you think ?

It looks good for right now. We can help you debug it when the rest of the wrappers are in place.

Thanks,

Matt

...

Thibault

Le dim. 13 déc. 2020 à 15:29, Matthew Knepley <[email protected]> a écrit :

...

On Sun, Dec 13, 2020 at 9:17 AM Mark Adams <[email protected]> wrote:

...

I don't think function pointers to PETSc (DS and DM) methods are not going to work in Fortran:

ierr = PetscDSAddBoundary(prob, DM_BC_NATURAL_RIEMANN, "inflow", "Face Sets", 0, 0, NULL, (void (*)(void)) PhysicsBoundary_Advect_Inflow, NULL, ALEN(inflowids), inflowids, phys);CHKERRQ(ierr);

You could write a funcs.c file that you call from your fortran code, like, call setBC1(prob,...,ierr)

and put PhysicsBoundary_Advect_Inflow and setBC1 in funcs.c, for instance.

Wrappers for functions that that function arguments are possible, but more involved. We have to create an internal struct that holds Fortran function pointers, and then call them with the right arguments. This has been done for SNESSetFunction() for instance, so we would start emulating that.

Thanks,

Matt

...

On Sun, Dec 13, 2020 at 5:32 AM Thibault Bridel-Bertomeu < [email protected]> wrote:

...

Good morning all,

Thank you Barry for your answer. I started adding some interfaces (PetscFVSetComponentName, PetscFVView & PetscFVSetType so far) and I have to say I think it is working quite well, but the prototypes of those functions are still quite "simple". I am stuck at how to implement the wrappers for PetscDSSetRiemannSolver and PetscDSSetContext though, especially on how to pass a function as an argument to PetscDSSetRiemannSolver ... Are there any similar functions that may already have their wrappers ?

Thank you very much,

Thibault

Le sam. 12 déc. 2020 à 23:28, Barry Smith <[email protected]> a écrit :

...

On Dec 12, 2020, at 2:59 PM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear Jed, dear Barry,

Thank you for the fast answers !

If I have any success I will make sure to make a pull request to provide some version of ex11 in Fortran.

Regarding the stubs, I admit I started looking in that direction to add the missing wrappers but I am not sure I fully understand the process yet. For each C function, I gotta provide a Fortran interface in a .h90 file as well as a C function that has a Fortran-like prototype and calls the C function - right ?

Yes,

However there are a few things I could not find / understand yet. For instance, it appears that for C functions that have character string arguments take an extra argument in their Fortran-like-prototype-wrapper, namely the length of the string. Is that passed automatically ? I couldn’t find where it could come from ...

This secret argument is put automatically by the Fortran compiler.

Another thing is for functions like PetscFVView. I guess the wrapping is less straightforward because I tried a quick something and it segfault’ed. I couldnt find the wrapper for DMView although there is such a routine in Fortran too. Could you please detail how to wrap such functions ?

PETSC_EXTERN void dmview_(DM *da,PetscViewer *vin,PetscErrorCode *ierr) { PetscViewer v; PetscPatchDefaultViewers_Fortran(vin,v); *ierr = DMView(*da,v); }

dm/interface/ftn-custom/zdmf.c

Thank you very much again,

Thibault Bridel-Bertomeu

Le sam. 12 déc. 2020 à 21:48, Barry Smith <[email protected]> a écrit :

...

PETSc Fortran interfaces are a combination of automatically generated and manually generated.

For any C PETSc function if the manual page begins with /*@ it generates the Fortran interface automatically (make allfortranstubs). If it begins /*@C then either the Fortran interface is done manually or is missing.

C functions that have character string arguments or function arguments (or a few other special cases) need to be manually provided. The automatically generated stubs go in the directory ftn-auto while manually generated ones go in the directory fin-custom.

Perhaps you could first generate a list of "missing" Fortran stubs and then for each stub determine why it is missing and if it can be provided. Some are likely easy to provide but a few (involving function arguments) will be more involved. Once you have all the stubs available translating ex11.c becomes straightforward.

Barry

On Dec 12, 2020, at 9:30 AM, Thibault Bridel-Bertomeu < [email protected]> wrote:

Dear all,

Is there somewhere a version of the TS tutorial ex11.c in Fortran ? I am looking into building in F90 (let's say that it is an unavoidable constraint) an unstructured 3D solver of the Euler equations using the "new" features of PETSc - mostly DMPlex & PetscFV - but I think there are some interfaces missing and I find it hard to find workarounds in Fortran. I would be grateful if anyone could please give me some pointers ...

Thank you very much in advance,

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

--

Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected]

-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>

-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>