Actual source code: test1.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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 rational function.\n\n";
13: #include <slepcfn.h>
15: int main(int argc,char **argv)
16: {
18: FN fn;
19: PetscInt i,na,nb;
20: PetscScalar x,y,yp,p[10],q[10],five=5.0,*pp,*qq;
21: char strx[50],str[50];
23: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
24: FNCreate(PETSC_COMM_WORLD,&fn);
26: /* polynomial p(x) */
27: na = 5;
28: p[0] = -3.1; p[1] = 1.1; p[2] = 1.0; p[3] = -2.0; p[4] = 3.5;
29: FNSetType(fn,FNRATIONAL);
30: FNRationalSetNumerator(fn,na,p);
31: FNView(fn,NULL);
32: x = 2.2;
33: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
34: FNEvaluateFunction(fn,x,&y);
35: FNEvaluateDerivative(fn,x,&yp);
36: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
37: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
38: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
39: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
41: /* inverse of polynomial 1/q(x) */
42: nb = 3;
43: q[0] = -3.1; q[1] = 1.1; q[2] = 1.0;
44: FNSetType(fn,FNRATIONAL);
45: FNRationalSetNumerator(fn,0,NULL); /* reset previous values */
46: FNRationalSetDenominator(fn,nb,q);
47: FNView(fn,NULL);
48: x = 2.2;
49: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
50: FNEvaluateFunction(fn,x,&y);
51: FNEvaluateDerivative(fn,x,&yp);
52: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
53: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
54: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
55: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
57: /* rational p(x)/q(x) */
58: na = 2; nb = 3;
59: p[0] = 1.1; p[1] = 1.1;
60: q[0] = 1.0; q[1] = -2.0; q[2] = 3.5;
61: FNSetType(fn,FNRATIONAL);
62: FNRationalSetNumerator(fn,na,p);
63: FNRationalSetDenominator(fn,nb,q);
64: FNSetScale(fn,1.2,0.5);
65: FNView(fn,NULL);
66: x = 2.2;
67: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
68: FNEvaluateFunction(fn,x,&y);
69: FNEvaluateDerivative(fn,x,&yp);
70: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
71: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
72: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
73: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
75: FNRationalGetNumerator(fn,&na,&pp);
76: FNRationalGetDenominator(fn,&nb,&qq);
77: PetscPrintf(PETSC_COMM_WORLD,"Coefficients:\n Numerator: ");
78: for (i=0;i<na;i++) {
79: SlepcSNPrintfScalar(str,50,pp[i],PETSC_FALSE);
80: PetscPrintf(PETSC_COMM_WORLD,"%s ",str);
81: }
82: PetscPrintf(PETSC_COMM_WORLD,"\n Denominator: ");
83: for (i=0;i<nb;i++) {
84: SlepcSNPrintfScalar(str,50,qq[i],PETSC_FALSE);
85: PetscPrintf(PETSC_COMM_WORLD,"%s ",str);
86: }
87: PetscPrintf(PETSC_COMM_WORLD,"\n");
88: PetscFree(pp);
89: PetscFree(qq);
91: /* constant */
92: FNSetType(fn,FNRATIONAL);
93: FNRationalSetNumerator(fn,1,&five);
94: FNRationalSetDenominator(fn,0,NULL); /* reset previous values */
95: FNView(fn,NULL);
96: x = 2.2;
97: SlepcSNPrintfScalar(strx,50,x,PETSC_FALSE);
98: FNEvaluateFunction(fn,x,&y);
99: FNEvaluateDerivative(fn,x,&yp);
100: SlepcSNPrintfScalar(str,50,y,PETSC_FALSE);
101: PetscPrintf(PETSC_COMM_WORLD," f(%s)=%s\n",strx,str);
102: SlepcSNPrintfScalar(str,50,yp,PETSC_FALSE);
103: PetscPrintf(PETSC_COMM_WORLD," f'(%s)=%s\n",strx,str);
105: FNDestroy(&fn);
106: SlepcFinalize();
107: return ierr;
108: }
110: /*TEST
112: test:
113: suffix: 1
114: nsize: 1
116: TEST*/