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