Actual source code: qeplin.c
slepc-3.8.2 2017-12-01
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2017, 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: */
10: /*
11: Various types of linearization for quadratic eigenvalue problem
12: */
14: #include <slepc/private/pepimpl.h>
15: #include "linearp.h"
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 ]
22: N1:
23: A = [ 0 I ] B = [ I 0 ]
24: [ -K -C ] [ 0 M ]
26: N2:
27: A = [ -K 0 ] B = [ C M ]
28: [ 0 I ] [ I 0 ]
30: S1:
31: A = [ 0 -K ] B = [-K 0 ]
32: [ -K -C ] [ 0 M ]
34: S2:
35: A = [ -K 0 ] B = [ C M ]
36: [ 0 M ] [ M 0 ]
38: H1:
39: A = [ K 0 ] B = [ 0 K ]
40: [ C K ] [-M 0 ]
42: H2:
43: A = [ 0 -K ] B = [ M C ]
44: [ M 0 ] [ 0 M ]
45: */
47: /* - - - N1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
49: PetscErrorCode MatCreateExplicit_Linear_N1A(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
50: {
52: PetscInt M,N,m,n,i,Istart,Iend;
53: Mat Id;
56: MatGetSize(ctx->M,&M,&N);
57: MatGetLocalSize(ctx->M,&m,&n);
58: MatCreate(PetscObjectComm((PetscObject)ctx->M),&Id);
59: MatSetSizes(Id,m,n,M,N);
60: MatSetFromOptions(Id);
61: MatSetUp(Id);
62: MatGetOwnershipRange(Id,&Istart,&Iend);
63: for (i=Istart;i<Iend;i++) {
64: MatSetValue(Id,i,i,1.0,INSERT_VALUES);
65: }
66: MatAssemblyBegin(Id,MAT_FINAL_ASSEMBLY);
67: MatAssemblyEnd(Id,MAT_FINAL_ASSEMBLY);
68: MatCreateTile(0.0,Id,1.0,Id,-ctx->dsfactor,ctx->K,-ctx->sfactor*ctx->dsfactor,ctx->C,A);
69: MatDestroy(&Id);
70: return(0);
71: }
73: PetscErrorCode MatCreateExplicit_Linear_N1B(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
74: {
76: PetscInt M,N,m,n,i,Istart,Iend;
77: Mat Id;
80: MatGetSize(ctx->M,&M,&N);
81: MatGetLocalSize(ctx->M,&m,&n);
82: MatCreate(PetscObjectComm((PetscObject)ctx->M),&Id);
83: MatSetSizes(Id,m,n,M,N);
84: MatSetFromOptions(Id);
85: MatSetUp(Id);
86: MatGetOwnershipRange(Id,&Istart,&Iend);
87: for (i=Istart;i<Iend;i++) {
88: MatSetValue(Id,i,i,1.0,INSERT_VALUES);
89: }
90: MatAssemblyBegin(Id,MAT_FINAL_ASSEMBLY);
91: MatAssemblyEnd(Id,MAT_FINAL_ASSEMBLY);
92: MatCreateTile(1.0,Id,0.0,Id,0.0,Id,ctx->sfactor*ctx->sfactor*ctx->dsfactor,ctx->M,B);
93: MatDestroy(&Id);
94: return(0);
95: }
97: /* - - - N2 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
99: PetscErrorCode MatCreateExplicit_Linear_N2A(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
100: {
102: PetscInt M,N,m,n,i,Istart,Iend;
103: Mat Id;
106: MatGetSize(ctx->M,&M,&N);
107: MatGetLocalSize(ctx->M,&m,&n);
108: MatCreate(PetscObjectComm((PetscObject)ctx->M),&Id);
109: MatSetSizes(Id,m,n,M,N);
110: MatSetFromOptions(Id);
111: MatSetUp(Id);
112: MatGetOwnershipRange(Id,&Istart,&Iend);
113: for (i=Istart;i<Iend;i++) {
114: MatSetValue(Id,i,i,1.0,INSERT_VALUES);
115: }
116: MatAssemblyBegin(Id,MAT_FINAL_ASSEMBLY);
117: MatAssemblyEnd(Id,MAT_FINAL_ASSEMBLY);
118: MatCreateTile(-1.0,ctx->K,0.0,Id,0.0,Id,1.0,Id,A);
119: MatDestroy(&Id);
120: return(0);
121: }
123: PetscErrorCode MatCreateExplicit_Linear_N2B(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
124: {
126: PetscInt M,N,m,n,i,Istart,Iend;
127: Mat Id;
130: MatGetSize(ctx->M,&M,&N);
131: MatGetLocalSize(ctx->M,&m,&n);
132: MatCreate(PetscObjectComm((PetscObject)ctx->M),&Id);
133: MatSetSizes(Id,m,n,M,N);
134: MatSetFromOptions(Id);
135: MatSetUp(Id);
136: MatGetOwnershipRange(Id,&Istart,&Iend);
137: for (i=Istart;i<Iend;i++) {
138: MatSetValue(Id,i,i,1.0,INSERT_VALUES);
139: }
140: MatAssemblyBegin(Id,MAT_FINAL_ASSEMBLY);
141: MatAssemblyEnd(Id,MAT_FINAL_ASSEMBLY);
142: MatCreateTile(ctx->sfactor,ctx->C,ctx->sfactor*ctx->sfactor,ctx->M,1.0,Id,0.0,Id,B);
143: MatDestroy(&Id);
144: return(0);
145: }
147: /* - - - S1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
149: PetscErrorCode MatCreateExplicit_Linear_S1A(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
150: {
154: MatCreateTile(0.0,ctx->K,-1.0,ctx->K,-1.0,ctx->K,-ctx->sfactor,ctx->C,A);
155: return(0);
156: }
158: PetscErrorCode MatCreateExplicit_Linear_S1B(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
159: {
163: MatCreateTile(-1.0,ctx->K,0.0,ctx->M,0.0,ctx->M,ctx->sfactor*ctx->sfactor,ctx->M,B);
164: return(0);
165: }
167: /* - - - S2 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
169: PetscErrorCode MatCreateExplicit_Linear_S2A(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
170: {
174: MatCreateTile(-1.0,ctx->K,0.0,ctx->M,0.0,ctx->M,ctx->sfactor*ctx->sfactor,ctx->M,A);
175: return(0);
176: }
178: PetscErrorCode MatCreateExplicit_Linear_S2B(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
179: {
183: MatCreateTile(ctx->sfactor,ctx->C,ctx->sfactor*ctx->sfactor,ctx->M,ctx->sfactor*ctx->sfactor,ctx->M,0.0,ctx->M,B);
184: return(0);
185: }
187: /* - - - H1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
189: PetscErrorCode MatCreateExplicit_Linear_H1A(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
190: {
194: MatCreateTile(1.0,ctx->K,0.0,ctx->K,ctx->sfactor,ctx->C,1.0,ctx->K,A);
195: return(0);
196: }
198: PetscErrorCode MatCreateExplicit_Linear_H1B(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
199: {
203: MatCreateTile(0.0,ctx->K,1.0,ctx->K,-ctx->sfactor*ctx->sfactor,ctx->M,0.0,ctx->K,B);
204: return(0);
205: }
207: /* - - - H2 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
209: PetscErrorCode MatCreateExplicit_Linear_H2A(MPI_Comm comm,PEP_LINEAR *ctx,Mat *A)
210: {
214: MatCreateTile(0.0,ctx->K,-1.0,ctx->K,ctx->sfactor*ctx->sfactor,ctx->M,0.0,ctx->K,A);
215: return(0);
216: }
218: PetscErrorCode MatCreateExplicit_Linear_H2B(MPI_Comm comm,PEP_LINEAR *ctx,Mat *B)
219: {
223: MatCreateTile(ctx->sfactor*ctx->sfactor,ctx->M,ctx->sfactor,ctx->C,0.0,ctx->C,ctx->sfactor*ctx->sfactor,ctx->M,B);
224: return(0);
225: }