Actual source code: test10.c

slepc-3.9.1 2018-05-02
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2018, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */

 11: static char help[] = "Test Phi functions.\n\n";

 13: #include <slepcfn.h>

 15: int main(int argc,char **argv)
 16: {
 18:   FN             phi0,phi1,phik,phicopy;
 19:   PetscInt       k;
 20:   PetscScalar    x,y,yp,tau,eta;
 21:   char           strx[50],str[50];

 23:   SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;

 25:   /* phi_0(x) = exp(x) */
 26:   FNCreate(PETSC_COMM_WORLD,&phi0);
 27:   FNSetType(phi0,FNPHI);
 28:   FNPhiSetIndex(phi0,0);
 29:   FNView(phi0,NULL);
 30:   x = 2.2;
 31:   SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
 32:   FNEvaluateFunction(phi0,x,&y);
 33:   FNEvaluateDerivative(phi0,x,&yp);
 34:   SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
 35:   PetscPrintf(PETSC_COMM_WORLD,"  f(%s)=%s\n",strx,str);
 36:   SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
 37:   PetscPrintf(PETSC_COMM_WORLD,"  f'(%s)=%s\n",strx,str);

 39:   /* phi_1(x) = (exp(x)-1)/x with scaling factors eta*phi_1(tau*x) */
 40:   FNCreate(PETSC_COMM_WORLD,&phi1);
 41:   FNSetType(phi1,FNPHI);  /* default index should be 1 */
 42:   tau = 0.2;
 43:   eta = 1.3;
 44:   FNSetScale(phi1,tau,eta);
 45:   FNView(phi1,NULL);
 46:   x = 2.2;
 47:   SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
 48:   FNEvaluateFunction(phi1,x,&y);
 49:   FNEvaluateDerivative(phi1,x,&yp);
 50:   SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
 51:   PetscPrintf(PETSC_COMM_WORLD,"  f(%s)=%s\n",strx,str);
 52:   SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
 53:   PetscPrintf(PETSC_COMM_WORLD,"  f'(%s)=%s\n",strx,str);

 55:   /* phi_k(x) with index set from command-line arguments */
 56:   FNCreate(PETSC_COMM_WORLD,&phik);
 57:   FNSetType(phik,FNPHI);
 58:   FNSetFromOptions(phik);

 60:   FNDuplicate(phik,PETSC_COMM_WORLD,&phicopy);
 61:   FNPhiGetIndex(phicopy,&k);
 62:   PetscPrintf(PETSC_COMM_WORLD,"Index of phi function is %D\n",k);
 63:   FNView(phicopy,NULL);
 64:   x = 2.2;
 65:   SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
 66:   FNEvaluateFunction(phicopy,x,&y);
 67:   FNEvaluateDerivative(phicopy,x,&yp);
 68:   SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
 69:   PetscPrintf(PETSC_COMM_WORLD,"  f(%s)=%s\n",strx,str);
 70:   SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
 71:   PetscPrintf(PETSC_COMM_WORLD,"  f'(%s)=%s\n",strx,str);

 73:   FNDestroy(&phi0);
 74:   FNDestroy(&phi1);
 75:   FNDestroy(&phik);
 76:   FNDestroy(&phicopy);
 77:   SlepcFinalize();
 78:   return ierr;
 79: }