1*6dd63270SBarry Smith #include <petscsnes.h>
2*6dd63270SBarry Smith #include <petsc/private/ftnimpl.h>
3*6dd63270SBarry Smith
4*6dd63270SBarry Smith #if defined(PETSC_HAVE_FORTRAN_CAPS)
5*6dd63270SBarry Smith #define snesgetconvergencehistory_ SNESGETCONVERGENCEHISTORY
6*6dd63270SBarry Smith #define snesrestoreconvergencehistory_ SNESRESTORECONVERGENCEHISTORY
7*6dd63270SBarry Smith #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
8*6dd63270SBarry Smith #define snesgetconvergencehistory_ snesgetconvergencehistory
9*6dd63270SBarry Smith #define snesrestoreconvergencehistory_ snesrestoreconvergencehistory
10*6dd63270SBarry Smith #endif
11*6dd63270SBarry Smith
snesgetconvergencehistory_(SNES * snes,F90Array1d * r,F90Array1d * fits,PetscInt * n,int * ierr PETSC_F90_2PTR_PROTO (ptrd1)PETSC_F90_2PTR_PROTO (ptrd2))12*6dd63270SBarry Smith PETSC_EXTERN void snesgetconvergencehistory_(SNES *snes, F90Array1d *r, F90Array1d *fits, PetscInt *n, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2))
13*6dd63270SBarry Smith {
14*6dd63270SBarry Smith PetscReal *hist;
15*6dd63270SBarry Smith PetscInt *its, N;
16*6dd63270SBarry Smith
17*6dd63270SBarry Smith CHKFORTRANNULLINTEGER(n);
18*6dd63270SBarry Smith *ierr = SNESGetConvergenceHistory(*snes, &hist, &its, &N);
19*6dd63270SBarry Smith if (*ierr) return;
20*6dd63270SBarry Smith *ierr = F90Array1dCreate(hist, MPIU_REAL, 1, N, r PETSC_F90_2PTR_PARAM(ptrd1));
21*6dd63270SBarry Smith if (*ierr) return;
22*6dd63270SBarry Smith *ierr = F90Array1dCreate(its, MPIU_INT, 1, N, fits PETSC_F90_2PTR_PARAM(ptrd2));
23*6dd63270SBarry Smith if (n) *n = N;
24*6dd63270SBarry Smith }
25*6dd63270SBarry Smith
snesrestoreconvergencehistory_(SNES * snes,F90Array1d * r,F90Array1d * fits,PetscInt * n,int * ierr PETSC_F90_2PTR_PROTO (ptrd1)PETSC_F90_2PTR_PROTO (ptrd2))26*6dd63270SBarry Smith PETSC_EXTERN void snesrestoreconvergencehistory_(SNES *snes, F90Array1d *r, F90Array1d *fits, PetscInt *n, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2))
27*6dd63270SBarry Smith {
28*6dd63270SBarry Smith *ierr = F90Array1dDestroy(r, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd1));
29*6dd63270SBarry Smith if (*ierr) return;
30*6dd63270SBarry Smith *ierr = F90Array1dDestroy(fits, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd2));
31*6dd63270SBarry Smith }
32