Actual source code: test1.c
slepc-3.5.2 2014-10-10
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: static char help[] = "Test rational function.\n\n";
24: #include <slepcfn.h>
28: int main(int argc,char **argv)
29: {
31: FN fn;
32: PetscInt na,nb;
33: PetscScalar x,y,yp,p[10],q[10],five=5.0;
34: char strx[50],str[50];
36: SlepcInitialize(&argc,&argv,(char*)0,help);
37: FNCreate(PETSC_COMM_WORLD,&fn);
39: /* polynomial p(x) */
40: na = 5;
41: p[0] = -3.1; p[1] = 1.1; p[2] = 1.0; p[3] = -2.0; p[4] = 3.5;
42: FNSetType(fn,FNRATIONAL);
43: FNSetParameters(fn,na,p,0,NULL);
44: FNView(fn,NULL);
45: x = 2.2;
46: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
47: FNEvaluateFunction(fn,x,&y);
48: FNEvaluateDerivative(fn,x,&yp);
49: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
50: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
51: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
52: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
54: /* inverse of polynomial 1/q(x) */
55: nb = 3;
56: q[0] = -3.1; q[1] = 1.1; q[2] = 1.0;
57: FNSetType(fn,FNRATIONAL);
58: FNSetParameters(fn,0,NULL,nb,q);
59: FNView(fn,NULL);
60: x = 2.2;
61: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
62: FNEvaluateFunction(fn,x,&y);
63: FNEvaluateDerivative(fn,x,&yp);
64: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
65: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
66: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
67: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
69: /* rational p(x)/q(x) */
70: na = 2; nb = 3;
71: p[0] = -3.1; p[1] = 1.1;
72: q[0] = 1.0; q[1] = -2.0; q[2] = 3.5;
73: FNSetType(fn,FNRATIONAL);
74: FNSetParameters(fn,na,p,nb,q);
75: FNView(fn,NULL);
76: x = 2.2;
77: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
78: FNEvaluateFunction(fn,x,&y);
79: FNEvaluateDerivative(fn,x,&yp);
80: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
81: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
82: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
83: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
85: /* constant */
86: FNSetType(fn,FNRATIONAL);
87: FNSetParameters(fn,1,&five,0,NULL);
88: FNView(fn,NULL);
89: x = 2.2;
90: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
91: FNEvaluateFunction(fn,x,&y);
92: FNEvaluateDerivative(fn,x,&yp);
93: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
94: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
95: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
96: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
98: FNDestroy(&fn);
99: SlepcFinalize();
100: return 0;
101: }