Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
6 : This file is part of SLEPc.
7 : SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9 : */
10 : /*
11 : Various types of linearization for quadratic eigenvalue problem
12 : */
13 :
14 : #include <slepc/private/pepimpl.h>
15 : #include "linear.h"
16 :
17 : /*
18 : Given the quadratic problem (l^2*M + l*C + K)*x = 0 the linearization is
19 : A*z = l*B*z for z = [ x ] and A,B defined as follows:
20 : [ l*x ]
21 :
22 : N:
23 : A = [-bK aI ] B = [ aI+bC bM ]
24 : [-aK -aC+bI ] [ bI aM ]
25 :
26 : S:
27 : A = [ bK aK ] B = [ aK-bC -bM ]
28 : [ aK aC-bM ] [-bM -aM ]
29 :
30 : H:
31 : A = [ aK -bK ] B = [ bM aK+bC ]
32 : [ aC+bM aK ] [-aM bM ]
33 : */
34 :
35 : /* - - - N - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
36 :
37 6 : PetscErrorCode MatCreateExplicit_Linear_NA(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
38 : {
39 6 : PetscInt M,N,m,n;
40 6 : Mat Id,T=NULL;
41 6 : PetscReal a=ctx->alpha,b=ctx->beta;
42 6 : PetscScalar scalt=1.0;
43 :
44 6 : PetscFunctionBegin;
45 6 : PetscCall(MatGetSize(ctx->M,&M,&N));
46 6 : PetscCall(MatGetLocalSize(ctx->M,&m,&n));
47 6 : PetscCall(MatCreateConstantDiagonal(PetscObjectComm((PetscObject)ctx->M),m,n,M,N,1.0,&Id));
48 6 : if (a!=0.0 && b!=0.0) {
49 0 : PetscCall(MatDuplicate(ctx->C,MAT_COPY_VALUES,&T));
50 0 : PetscCall(MatScale(T,-a*ctx->dsfactor*ctx->sfactor));
51 0 : PetscCall(MatShift(T,b));
52 : } else {
53 6 : if (a==0.0) { T = Id; scalt = b; }
54 5 : else { T = ctx->C; scalt = -a*ctx->dsfactor*ctx->sfactor; }
55 : }
56 6 : PetscCall(MatCreateTile(-b*ctx->dsfactor,ctx->K,a,Id,-ctx->dsfactor*a,ctx->K,scalt,T,A));
57 6 : PetscCall(MatDestroy(&Id));
58 6 : if (a!=0.0 && b!=0.0) PetscCall(MatDestroy(&T));
59 6 : PetscFunctionReturn(PETSC_SUCCESS);
60 : }
61 :
62 6 : PetscErrorCode MatCreateExplicit_Linear_NB(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
63 : {
64 6 : PetscInt M,N,m,n;
65 6 : Mat Id,T=NULL;
66 6 : PetscReal a=ctx->alpha,b=ctx->beta;
67 6 : PetscScalar scalt=1.0;
68 :
69 6 : PetscFunctionBegin;
70 6 : PetscCall(MatGetSize(ctx->M,&M,&N));
71 6 : PetscCall(MatGetLocalSize(ctx->M,&m,&n));
72 6 : PetscCall(MatCreateConstantDiagonal(PetscObjectComm((PetscObject)ctx->M),m,n,M,N,1.0,&Id));
73 6 : if (a!=0.0 && b!=0.0) {
74 0 : PetscCall(MatDuplicate(ctx->C,MAT_COPY_VALUES,&T));
75 0 : PetscCall(MatScale(T,b*ctx->dsfactor*ctx->sfactor));
76 0 : PetscCall(MatShift(T,a));
77 : } else {
78 6 : if (b==0.0) { T = Id; scalt = a; }
79 1 : else { T = ctx->C; scalt = b*ctx->dsfactor*ctx->sfactor; }
80 : }
81 6 : PetscCall(MatCreateTile(scalt,T,b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,b,Id,a*ctx->sfactor*ctx->sfactor*ctx->dsfactor,ctx->M,B));
82 6 : PetscCall(MatDestroy(&Id));
83 6 : if (a!=0.0 && b!=0.0) PetscCall(MatDestroy(&T));
84 6 : PetscFunctionReturn(PETSC_SUCCESS);
85 : }
86 :
87 : /* - - - S - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
88 :
89 3 : PetscErrorCode MatCreateExplicit_Linear_SA(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
90 : {
91 3 : Mat T=NULL;
92 3 : PetscScalar scalt=1.0;
93 3 : PetscReal a=ctx->alpha,b=ctx->beta;
94 :
95 3 : PetscFunctionBegin;
96 3 : if (a!=0.0 && b!=0.0) {
97 0 : PetscCall(MatDuplicate(ctx->C,MAT_COPY_VALUES,&T));
98 0 : PetscCall(MatScale(T,a*ctx->dsfactor*ctx->sfactor));
99 0 : PetscCall(MatAXPY(T,-b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,UNKNOWN_NONZERO_PATTERN));
100 : } else {
101 3 : if (a==0.0) { T = ctx->M; scalt = -b*ctx->dsfactor*ctx->sfactor*ctx->sfactor; }
102 2 : else { T = ctx->C; scalt = a*ctx->dsfactor*ctx->sfactor; }
103 : }
104 3 : PetscCall(MatCreateTile(b*ctx->dsfactor,ctx->K,a*ctx->dsfactor,ctx->K,a*ctx->dsfactor,ctx->K,scalt,T,A));
105 3 : if (a!=0.0 && b!=0.0) PetscCall(MatDestroy(&T));
106 3 : PetscFunctionReturn(PETSC_SUCCESS);
107 : }
108 :
109 3 : PetscErrorCode MatCreateExplicit_Linear_SB(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
110 : {
111 3 : Mat T=NULL;
112 3 : PetscScalar scalt=1.0;
113 3 : PetscReal a=ctx->alpha,b=ctx->beta;
114 :
115 3 : PetscFunctionBegin;
116 3 : if (a!=0.0 && b!=0.0) {
117 0 : PetscCall(MatDuplicate(ctx->C,MAT_COPY_VALUES,&T));
118 0 : PetscCall(MatScale(T,-b*ctx->dsfactor*ctx->sfactor));
119 0 : PetscCall(MatAXPY(T,a*ctx->dsfactor,ctx->K,UNKNOWN_NONZERO_PATTERN));
120 : } else {
121 3 : if (b==0.0) { T = ctx->K; scalt = a*ctx->dsfactor; }
122 1 : else { T = ctx->C; scalt = -b*ctx->dsfactor*ctx->sfactor; }
123 : }
124 3 : PetscCall(MatCreateTile(scalt,T,-b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,-b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,-a*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,B));
125 3 : if (a!=0.0 && b!=0.0) PetscCall(MatDestroy(&T));
126 3 : PetscFunctionReturn(PETSC_SUCCESS);
127 : }
128 :
129 : /* - - - H - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
130 :
131 2 : PetscErrorCode MatCreateExplicit_Linear_HA(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
132 : {
133 2 : Mat T=NULL;
134 2 : PetscScalar scalt=1.0;
135 2 : PetscReal a=ctx->alpha,b=ctx->beta;
136 :
137 2 : PetscFunctionBegin;
138 2 : if (a!=0.0 && b!=0.0) {
139 0 : PetscCall(MatDuplicate(ctx->C,MAT_COPY_VALUES,&T));
140 0 : PetscCall(MatScale(T,a*ctx->dsfactor*ctx->sfactor));
141 0 : PetscCall(MatAXPY(T,b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,UNKNOWN_NONZERO_PATTERN));
142 : } else {
143 2 : if (a==0.0) { T = ctx->M; scalt = b*ctx->dsfactor*ctx->sfactor*ctx->sfactor; }
144 1 : else { T = ctx->C; scalt = a*ctx->dsfactor*ctx->sfactor; }
145 : }
146 2 : PetscCall(MatCreateTile(a*ctx->dsfactor,ctx->K,-b*ctx->dsfactor,ctx->K,scalt,T,a*ctx->dsfactor,ctx->K,A));
147 2 : if (a!=0.0 && b!=0.0) PetscCall(MatDestroy(&T));
148 2 : PetscFunctionReturn(PETSC_SUCCESS);
149 : }
150 :
151 2 : PetscErrorCode MatCreateExplicit_Linear_HB(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
152 : {
153 2 : Mat T=NULL;
154 2 : PetscScalar scalt=1.0;
155 2 : PetscReal a=ctx->alpha,b=ctx->beta;
156 :
157 2 : PetscFunctionBegin;
158 2 : if (a!=0.0 && b!=0.0) {
159 0 : PetscCall(MatDuplicate(ctx->C,MAT_COPY_VALUES,&T));
160 0 : PetscCall(MatScale(T,b*ctx->dsfactor*ctx->sfactor));
161 0 : PetscCall(MatAXPY(T,a*ctx->dsfactor,ctx->K,UNKNOWN_NONZERO_PATTERN));
162 : } else {
163 2 : if (b==0.0) { T = ctx->K; scalt = a*ctx->dsfactor; }
164 1 : else { T = ctx->C; scalt = b*ctx->dsfactor*ctx->sfactor; }
165 : }
166 2 : PetscCall(MatCreateTile(b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,scalt,T,-a*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,b*ctx->dsfactor*ctx->sfactor*ctx->sfactor,ctx->M,B));
167 2 : if (a!=0.0 && b!=0.0) PetscCall(MatDestroy(&T));
168 2 : PetscFunctionReturn(PETSC_SUCCESS);
169 : }
|