Sistemas Acoplados Edp PDF

UNIVERSIDAD POLITÉCNICA DE VALENCIA
DEPARTAMENTO DE MATEMÁTICA APLICADA
SOLUCIONES NUMÉRICAS
ESTABLES DE SISTEMAS
ACOPLADOS MIXTOS DE
ECUACIONES EN DERIVADAS
PARCIALES
TESIS DOCTORAL
Presentada por:
Marı́a Consuelo Casabán Bartual
Dirigida por:
Dr.D. Lucas Jódar Sánchez
Valencia, abril de 2002

D. LUCAS JÓDAR SÁNCHEZ, Catedrático de Universidad del Departamento de Ma-
temática Aplicada de la Universidad Politécnica de Valencia
CERTIFICA
Que la presente memoria Soluciones Numéricas Estables de Sistemas Acoplados Mixtos de
Ecuaciones en Derivadas Parciales ha sido realizada bajo su dirección, en el Departamen-
to de Matemática Aplicada de la Universidad Politécnica de Valencia, por la licenciada
MARÍA CONSUELO CASABÁN BARTUAL y constituye su tesis para optar al grado de
Doctor en Ciencias Matemáticas.
Para que conste, en cumplimiento de la legislación vigente, autoriza la presentación de
la referida tesis doctoral ante la comisión de Doctorado de la Universidad Politécnica de
Valencia, firmando el presente certificado.
Valencia, abril de 2002
Fdo: Lucas Jódar Sánchez

A la meua famı́lia, en especial a ma mare
per haver-me animat a realitzar els estudis
de 3er cicle
“ Tended a ser un poco aprendices de
todo, para vuestro bien, y maestros en
algo, para bien de los demás ”.
“ No existe problema alguno en ingenierı́a

que no esté subordinado a una imagen o
representación de carácter matemático,
cuya solución, fácil o difı́cil, significará el
total dominio del campo correspondiente.
La Matemática es la ciencia de los esque-
mas aplicables a tales problemas ”.
P. Puig Adam
Agradecimientos
Cuando hace unos años empecé mi andadura por el difı́cil camino de la investigación, siendo
becaria de este departamento y dirigida por mi apreciado profesor Dr. Lucas Jódar, nunca
imaginé, que tal dı́a como hoy, llegarı́a a escribir estas letras, que de algún modo concluyen
mis estudios de tercer ciclo y que serán el principio de una nueva etapa.
No sin antes terminar, quisiera hacer constar en estas lı́neas, mi más sincero agradecimiento
al Dr. Lucas Jódar, mi director de tesis, por su ayuda y constancia en la realización de la
misma, y por transmitirme sus innumerables conocimientos sobre cómo investigar y ser un
buen docente.
Ası́ también, agradezco a mis amigos Dolors, Lola y Emilio su constante apoyo y ánimo.
i
Resumen
Esta memoria trata sobre la construcción de soluciones numéricas estables de sistemas

parabólicos e hiperbólicos acoplados. Las etapas caracterı́sticas de esta memoria son: la
construcción de soluciones discretas utilizando diferencias finitas y una técnica de sepa-
ración de variables discreta, el estudio de la estabilidad y la consistencia de la solución
calculada, y el empleo de un método de proyecciones para extender los resultados obteni-
dos a una clase más general de funciones de valores iniciales.
Mediante la aplicación de un método de separación de variables discreto, la solución numéri-

ca propuesta a los problemas, es la solución exacta de un sistema en diferencias acoplado,
que se obtiene de la discretización en diferencias finitas del sistema acoplado en derivadas
parciales continuo.
Las condiciones de contorno de los problemas aquı́ tratados son acopladas y de tipo no-
Dirichlet.
Nuestro enfoque metodológico es alternativo frente al tratamiento algebraico más tradicio-

nal que escribe el esquema matricialmente, y ofrece la ventaja de no tener que resolver los
sistemas algebraicos de gran tamaño con bloques matriciales que aparecen en el método
de diferencias finitas estándar, gracias al empleo de un método de separación de variables
discreto.
Las técnicas de desacoplamiento no son aplicables salvo para hipótesis excesivamente res-
trictivas. Aún desacoplando la ecuación en derivadas parciales, las condiciones de contorno
no tienen porqué estarlo, por lo que las técnicas de desacoplamiento son poco convenientes.
Los problemas tratados modelizan, entre otros, problemas de difusión, conducción nerviosa
y problemas del armamento (capı́tulo 2), calentamiento por microondas, óptica, cardiologı́a
y flujos del suelo (capı́tulo 3).
iii
Resum
Aquesta memòria tracta sobre la construcció de solucions numèriques estables de sistemes

parabòlics e hiperbòlics acoblats. Les etapes caracterı́stiques d’aquesta memòria són: la
construcció de solucions discretes utilitzant diferències finites i una tècnica de separació de
variables discreta, l’estudi de l’estabilitat i la consistència de la solució calculada, i l’ ús
d’un mètode de projeccions per a estendre els resultats obtinguts a una classe més general
de funcions de valors inicials.
Mitjançant l’aplicació d’un mètode de separació de variables discret, la solució numèrica

proposta als problemes, és la solució exacta d’un sistema en diferències acoblat, que s’obté de
la discretització en diferències finites del sistema acoblat en derivades parcials continu.
Les condicions de contorn dels problemes acı́ tractats són acoblades i de tipus no-Dirichlet.
El nostre enfocament metodològic és alternatiu front al tractament algebraic més tradi-
cional que escriu l’esquema matricialment, i ofereix l’avantatge de no haver de resoldre
els sistemes algebraics de gran tamany amb blocs matricials que apareixen en el mètode
de diferències finites estàndard, gràcies a l’ ús d’un mètode de separació de variables discret.
Les tècniques de desacoblament no són aplicables excepte per a hipòtesis excessivament res-
trictives. Fins i tot desacoblant l’equació en derivades parcials, les condicions de contorn no
tenen per què estar-ho, per la qual cosa les tècniques de desacoblament són poc convenients.
Els problemes tractats modelen, entre altres, problemes de difusió, conducció nerviosa i
problemes de l’armament (capı́tol 2), calfament per microones, òptica, cardiologia i fluixes
del sòl (capı́tol 3).
v
Abstract
This memory deals with the construction of stable numerical solutions of coupled parabolic
and hyperbolic systems. The characteristic stages of this memory are: the construction of
discrete solutions using difference schemes and a discrete separation of variables method,
the study of the stability and the consistency of the calculated solution, and the use of a
method of projections in order to extend the obtained results to a more general class of
initial value functions.
By the application of a discrete separation of variables method, the proposed numerical

solution to the problems, is the exact solution of a coupled difference system, which is ob-
tained from the discretization in finite differences of the coupled mixed partial differential
system.
The boundary conditions of the problems which have been dealt here are coupled and of
non-Dirichlet type.
Our methodological approach is alternative to the more traditional algebraic treatment

which writes the scheme in matrix form, and offers the advantage of not having to solve
the big-sized algebraic systems with matrix blocks appearing in the standard difference
method, thanks to the use of a discrete separation of variables method.
The uncoupled techniques aren’t applied except for excessively restrictive hypotheses. Even
uncoupling the partial differential equation, the boundary conditions don’t necessarily have
to be so, hence the uncoupled techniques aren’t very suitable.
The treated problems model, among others, diffusion problems, nerve conduction and ar-
mament models (chapter 2), microwave heating processes, optics, cardiology and soil flows
(chapter 3).
vii
Índice general
1. Introducción, motivación y preliminares 1

1.1. Introducción y motivación . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2. Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
2. Sistemas parabólicos 17
2.1. Caso particular . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.1.1. Introducción . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.1.2. La discretización del problema . . . . . . . . . . . . . . . . . . . . . 18
2.1.3. El problema de contorno discreto . . . . . . . . . . . . . . . . . . . . 19
2.1.4. Consistencia y estabilidad de las soluciones . . . . . . . . . . . . . . 26
2.1.5. El problema mixto discreto . . . . . . . . . . . . . . . . . . . . . . . 31
2.1.6. El método de las proyecciones . . . . . . . . . . . . . . . . . . . . . . 40
2.1.7. Ejemplos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
2.2. Caso general . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
2.2.1. Introducción . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
2.2.2. El problema de contorno discreto . . . . . . . . . . . . . . . . . . . . 62
2.2.3. El problema mixto discreto . . . . . . . . . . . . . . . . . . . . . . . 72
2.2.4. Ejemplos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78
3. Sistemas hiperbólicos 85
3.1. Introducción . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
3.2. La discretización del problema . . . . . . . . . . . . . . . . . . . . . . . . . 86
3.3. El problema de contorno discreto . . . . . . . . . . . . . . . . . . . . . . . . 86
3.4. El problema mixto discreto . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
3.5. Consistencia y estabilidad de las soluciones . . . . . . . . . . . . . . . . . . 98
3.6. El método de las proyecciones . . . . . . . . . . . . . . . . . . . . . . . . . . 101
3.7. Ejemplos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105
Bibliografı́a 109
ix
Capı́tulo 1
Introducción, motivación y
preliminares
0.4pt0.4pt 0pt0.4pt
1.1. Introducción y motivación

Esta memoria trata sobre la construción de soluciones numéricas estables de sistemas pa-
rabólicos e hiperbólicos acoplados con condiciones de contorno acopladas. Las etapas carac-
terı́sticas de esta memoria son: la construcción de soluciones discretas utilizando diferencias
finitas y una técnica de separación de variables discreta, el estudio de la estabilidad y la
consistencia de la solución calculada, y el empleo de un método de proyecciones para ex-
tender los resultados obtenidos a una clase más general de funciones de valores iniciales.
En este sentido, la memoria es una continuación de las tesis [28, 41], aunque en este trabajo
la ecuación matricial en derivadas parciales de los sistemas parabólicos considerados es más
general y además se incluye el estudio de sistemas hiperbólicos.
Debido a la complejidad que involucra el acoplamiento de los sistemas matriciales, los

métodos usuales para la resolución de estos problemas transforman el sistema acoplado en
un nuevo sistema desacoplado para poder abordarlo escalarmente [9, 42]. Estas técnicas
de desacoplamiento, a pesar de ser tan utilizadas, tienen desventajas bien conocidas tales
como: asumir hipótesis innecesarias, incrementar el orden de derivabilidad del sistema des-
acoplado, ser técnicas de aplicación limitada ya que la matriz de coeficientes tiene que ser
simétrica, y otras [8]. Además el desacoplamiento de la ecuación en derivadas parciales no
permite desacoplar las condiciones de contorno.
Nuestra oferta metodológica es alternativa frente al tratamiento tradicional que trata de
resolver el problema algebraico discretizado involucrando matrices por bloques cuyo espec-
tro es muy difı́cil de conocer o controlar, y ofrece la ventaja de no tener que resolver los
sistemas algebraicos de gran tamaño con bloques matriciales que aparecen en el método
de diferencias finitas estándar, gracias al empleo de un método de separación de variables
discreto. Además la construcción de las soluciones de los problemas planteados se hace de
forma exacta.
1
2 1.1 Introducción y motivación
La memoria está estructurada en dos capı́tulos. En el capı́tulo 2 se construyen soluciones

numéricas estables para sistemas parabólicos, considerando primeramente, condiciones de
contorno acopladas particulares y después de tipo más general. El estudio de estos dos casos
se ha desarrollado en los apartados 2.1 y 2.2. Los sistemas parabólicos estudiados son del
tipo:

ut (x, t) − Auxx (x, t) − Bu(x, t) = 0, 0 < x < 1, t > 0, 

u(0, t) = 0, t > 0,

, (1.1)
Cu(1, t) + Dux (1, t) = 0, t > 0, 

u(x, 0) = F (x), 0 ≤ x ≤ 1,


ut (x, t) − Auxx (x, t) − Bu(x, t) = 0, 0 < x < 1, t > 0, 

A1 u(0, t) + B1 ux (0, t) = 0, t > 0,

, (1.2)
A2 u(1, t) + B2 ux (1, t) = 0, t > 0, 

u(x, 0) = F (x), 0 ≤ x ≤ 1,

siendo A, B, C, D, Ak , Bk , k = 1, 2, matrices cuadradas en Cs×s , y la incógnita u y F

vectores en Cs .
Los sistemas mixtos difusivos de ecuaciones en derivadas parciales acoplados del tipo (1.1)
y (1.2) aparecen en el estudio de muchos campos de la ciencia tales como, geomecánica
[42]; fı́sica-quı́mica; mecánica [36]; en problemas de dispersión de mecánica cuántica [2, 31];
en la modelización termoelastoplástica; en problemas de conducción nerviosa [30, 17, 30];
en modelos del armamento [14]; en el encendido de una única componente de un gas no
reactivo en un recipiente cilı́ndrico cerrado con conservación de masa [27]; o en el estudio
de la muerte cardı́aca repentina como consecuencia de la fibrilación ventricular [48].
Las soluciones numéricas discretas para la ecuación
ut (x, t) − Auxx (x, t) = 0 , (1.3)
considerándose las mismas condiciones de contorno que en los problemas (1.1) y (1.2), han
sido tratadas en [26], y sus soluciones analı́tico-numéricas en [23].
La enorme profusión de los procesos de difusión se debe en parte, a la gran sencillez de las
dos ecuaciones que se requieren para poder establecer la ecuación de difusión. La primera
de ellas, denominada Primera ley de difusión, afirma que“la materia tiende a repartirse
de forma natural en el interior del sistema que lo contiene”. Aquı́ el concepto de materia
debe entenderse en su sentido más amplio, ya que nos referiremos por igual a las moléculas
de una sustancia, que tienden a desplazarse de las zonas de mayor concentración hacia las
zonas de menor concentración; o a las cargas, que se desplazan en contra de un gradiente
de potencial eléctrico; o al calor, que se reparte en un sólido.
La segunda ecuación de difusión, denominada Ley de la conservación de la materia/energı́a,
o también Ecuación de continuidad, es una ecuación de cumplimiento más general que la
Capı́tulo 1. Introducción, motivación y preliminares 3
anterior pues afirma que, “dado un volumen arbitrario, la variación de la cantidad de ma-
teria es igual a lo que entra en el volumen, menos lo que sale, más lo que se crea en su
interior, menos lo que se destruye”.
Como consecuencia de lo que se ha comentado, se concluye que un problema de difusión es
de naturaleza intrı́nsecamente escalar, debido a que no concierne más que a la descripción
de una determinada sustancia cuya caracterización se realiza por medio de una magni-
tud escalar como puede ser la concentración. A pesar de que no existe ninguna magnitud
conocida que verifique la ecuación de difusión, sı́ que existe relación entre la difusión de
distintas magnitudes en un mismo medio. Esto nos va a permitir expresar el modelo de
un problema mediante una ecuación de difusión matricial, donde cada entrada del vector
incógnita será una de las magnitudes cuya difusión simultánea estamos estudiando.
El acoplamiento entre varios problemas escalares se puede expresar de forma más elegante
mediante una formulación matricial. Atendiendo a [42, 34, 46] podemos modelar diferentes
fenómenos de conducción mediante el siguiente sistema general, transcrito en su versión
unidimensional:

∂ ∂ ∂ ∂U
(AU ) + (BU ) − D + EU + F = 0 , (1.4)
∂t ∂x ∂x ∂x
x ∈ Ω, t > 0.
El vector incógnita U es el vector de variables de estado, que en el caso más general con-
tendrá la temperatura, el potencial eléctrico, etc. A, B, D y E son matrices de coeficientes
y F es un vector de términos fuente.
Notemos que si A es invertible, D no depende de x, B = 0 y F = 0, entonces la ecuación
de difusión matricial (1.4) se transforma en la ecuación de los problemas (1.1) y (1.2) que
se estudian en el capı́tulo 2 de esta memoria, pero particularizada a la matriz A−1 D. Estas
condiciones corresponden a problemas fı́sicos en los que no tenemos convección (B = 0), no
hay términos fuente (F = 0), existen especies que se transforman en otras especies (E 6= 0),
y el medio es homogéneo (D es constante). Estas condiciones son las que se estudian en la
presente memoria.
Esta memoria concluye con el capı́tulo 3 en el que construyen soluciones numéricas estables
de sistemas hiperbólicos acoplados del tipo:

Auxx (x, t) − utt (x, t) = 0, 0 < x < 1, t > 0, 

u(0, t) = 0, t > 0,



Bu(1, t) + Cux (1, t) = 0, t > 0, , (1.5)
u(x, 0) = F (x), 0 ≤ x ≤ 1,




ut (x, 0) = V (x), 0 ≤ x ≤ 1.

siendo A, B, C matrices cuadradas en Cs×s , y la incógnita u y F , V vectores en Cs .

Los sistemas hiperbólicos fuertemente acoplados de ecuaciones en derivadas parciales del
tipo (1.5) surgen en procesos de calentamiento por microondas [32], [15], óptica [10], car-
diologı́a [47], flujos del suelo [42, 49], y otros.
0.4pt0.4pt 0pt0.4pt
4 1.2 Preliminares
1.2. Preliminares
A lo largo del trabajo vamos a utilizar una serie de definiciones, notación y resultados que
creemos conveniente enunciar en este apartado de preliminares para facilitar su lectura.
1/2
Definición 1.1 Si v es un vector en Cm , denotamos por kvk2 = v H v la norma
euclı́dea usual de v. Si A es una matriz en Cn×m , su norma 2 la denotamos por kAk,
y está definida por [13, pág. 56]
kA vk2
kAk = sup .
v6=0 kvk2
Definición 1.2 Si v es un vector en Cm , su norma infinito la denotamos por kvk∞ , y

está definida por [13, pág. 52]
kvk∞ = máx |vj | .
1≤j≤m
Definición 1.3 Una matriz A se dice convergente si la sucesión {An } tiende a la matriz
nula cuando n −→ ∞.
El conjunto de todos los valores propios de una matriz A en Cm×m lo denotaremos por
σ(A), y el radio espectral de A definido por el conjunto
máx {|z| ; z ∈ σ(A)} ,
lo denotaremos por ρ(A).

El núcleo de una matriz A, se denota por Nuc A, y la imagen de una matriz A por Im A.
Una caracterización del concepto de matriz convergente nos la proporciona el lema 5.2 de
[4, pág. 163]:
LEMA 1.1
A es convergente ⇐⇒ ρ(A) < 1 .
Definición 1.4 [3]. Dadas dos sucesiones {an } y {bn } tales que bn ≥ 0 para todo n, escri-
biremos
an = O(bn ) ,
si existe una constante M > 0 tal que |an | ≤ M bn para todo n. Por tanto, una ecuación de
la forma
an = cn + O(bn )
significa que
an − cn = O(bn ) .
En este trabajo utilizaremos un tipo de inversa generalizada, la inversa generalizada de

Moore-Penrose, que nos interesará para expresar en forma cerrada, finita y computable,
tanto la condición de compatibilidad como la solución general del sistema algebraico:
Ax = b, A ∈ Cm×n , x ∈ Cn , b ∈ Cm . (1.6)
Definición 1.5 Si A ∈ Cm×n , entonces su inversa generalizada Moore-Penrose se denota

por A† y se define como la única matriz que cumple:
(i) AA† A = A
(ii) A† AA† = A†
H
(iii) AA† = AA†
H
(iv) A† A = A† A
donde AH denota la traspuesta conjugada de A.
Un estudio más detallado de las propiedades y aplicaciones de este concepto se puede ver
en [6] y [39].
Para calcular las soluciones de (1.6) necesitaremos métodos eficaces para el cálculo de A† ,
como por ejemplo el programa Matlab, véase [16].
TEOREMA 1.1 (Teorema de Mitra) [39]. Sea A ∈ Cm×n y sea AG ∈ Cn×m una
inversa generalizada de A. Sea b ∈ Cm , entonces el sistema
Ax = b ,
es compatible (i.e. tiene solución) si y sólo si
AAG b = b .
Además en este caso, la solución general de Ax = b viene dada por
x = AG b + (I − AG A)z ,
donde z es un vector arbitrario de Cn .
Observación 1.1 [6]. Los siguientes resultados relacionan el núcleo y la imagen de una
matriz:
Nuc B = Im I − B † B ; Im B = Nuc I − BB † .
Definición 1.6 [13, pág. 311]. Un subespacio E de Cm es invariante por la matriz A ∈

Cm×m si
∀x ∈ E =⇒ Ax ∈ E .
Observación 1.2 Si tenemos v ∈ Cm tal que v ∈ Nuc B y además Nuc B es un subespacio

invariante por A entonces Aj v ∈ Nuc B, ∀j ≥ 0, i.e. BAj v = 0, ∀j ≥ 0. En efecto, puesto
que v ∈ Nuc B y Nuc B es un subespacio invariante por A tenemos que Av ∈ Nuc B.
Si ahora utilizamos de nuevo que Nuc B es un subespacio invariante por A pero tomando
Av ∈ Nuc B se obtiene que A2 v ∈ Nuc B. Por recurrencia, razonando análogamente se
obtiene que Aj v ∈ Nuc B, ∀j ≥ 0.
6 1.2 Preliminares
Observación 1.3 Como consecuencia de la observación 1.1 se obtienen las siguientes equi-
valencias:

Nuc B subespacio invariante por la matriz A ⇐⇒ BA I − B † B = 0 ,

Im B subespacio invariante por la matriz A ⇐⇒ I − BB † AB = 0 .
LEMA 1.2 [23]. Sean M y N matrices en Cs×s , entonces:
Nuc M ∩ Nuc N
h i† h i
† † †
= Im I −M M I − N I −M M N I −M M .
Demostración.
⊂
Si v ∈ Nuc M ∩ Nuc N entonces M v = 0 y N v =0. Por el teorema 2.3.2 de [1, pág. 24]
† †
obtenemos que M v = 0 implica v ∈ Im I − M M , i.e., v = I − M M d, donde d es un
vector arbitrario en Cs . Multiplicando v por N obtenemos:

0 = N v = N I − M †M d ,
i.e., d ∈ Nuc N I − M † M . Aplicando nuevamente el teorema 2.3.2 obtenemos:

h i† h i
† †
d= I − N I −M M N I −M M c, c ∈ Cs .
Multiplicando d por I − M † M :

v = I − M †M d
h i† h i
† † †
= I −M M I − N I −M M N I −M M c,
i.e., h i† h i

† † †
v ∈ Im I −M M I − N I −M M N I −M M .
Por tanto:
Nuc M ∩ Nuc N
h i† h i
† † †
⊂ Im I −M M I − N I −M M N I −M M .
⊃ n n † oo
Sea v ∈ Im I − M † M I − N I − M † M N I − M †M

.
s
Entonces para algún z ∈ C se tiene:
h i† h i
† † †
v = I −M M I − N I −M M N I −M M z.
De aquı́ multiplicando por M y N respectivamente y utilizando
M = M M † M,
por ser M † la inversa de Moore-Penrose, se sigue:

n † o
M v = M − M M †M I − N I − M †M N I − M †M

z = 0,
n † o
N v = N I − M †M − N I − M †M N I − M †M N I − M †M

z = 0 . De esta forma
v ∈ Nuc M ∩ Nuc N .
√
En todo el trabajo denotaremos por i a la unidad imaginaria, i = −1.
Definición 1.7 Toda matriz compleja A en Cm×m se puede descomponer de la siguiente

forma
A = A1 + i A 2 ,
donde A1 y A2 representan a la parte real e imaginaria de A respectivamente, las cuales se
definen por
A + AH A − AH
A1 = ; A2 = ,
2 2i
donde AH denota la traspuesta conjunda de A.
Las matrices A1 y A2 son matrices hermı́ticas, es decir:
AH
k = Ak , para k = 1, 2.
TEOREMA 1.2 [4]. Para una matriz A arbitraria se verifica

q
kAk = ρ (AH A) .
COROLARIO 1.1 [4]. Si P es una matriz hermı́tica (o real, simétrica) entonces
kP k = ρ(P ) .
Observación 1.4 [4]. Para una matriz arbitraria se verifica
ρ(A) ≤ kAk .
Observación 1.5 Sea P una matriz diagonalizable y Q una matriz invertible tal que
Q−1 P Q es una matriz diagonal, entonces kP k ≤ Q−1 kQkρ(P ).

TEOREMA 1.3 (Teorema de Bromwich) [33, pág. 389]. Sea A cualquier matriz com-
pleja. Sean µ y M el menor y el mayor valor propio de la matriz hermı́tica 12A + AH ,

1
y ν y N el menor y el mayor valor propio de la matriz hermı́tica 2i A − AH . Si λ es
cualquier valor propio de A, entonces,
µ ≤ Re(λ) ≤ M , ν ≤ Im(λ) ≤ N .
8 1.2 Preliminares
TEOREMA 1.4 [44, pág. 25]. Sean A y B matrices hermı́ticas, n × n, y A e = A + B.

Sean λ1 ≥ λ2 ≥ . . . ≥ λn los valores propios de A, y λ1 ≥ λ2 ≥ . . . ≥ λn los valores propios
e e e
de A.
e Si µn es el valor propio de B más pequeño, entonces
ei ≥ λi + µn ,
λ i = 1, . . . , n .
Si A es una matriz hermı́tica denotaremos
λmı́n (A) = mı́n (σ(A)) y λmáx (A) = máx (σ(A)) .
TEOREMA 1.5 (Teorema de Gershgorin) [4, pág. 127]. Dada una matriz A = [aij ]
se verifica que σ(A) está contenido en la unión de los discos
 
 X 
Fi = z ∈ C : |z − aii | ≤ |aij | , 1 ≤ i ≤ n ,
 
j6=i
y en la unión de los discos

 
 X 
Gi = z ∈ C : |z − aii | ≤ |aji | , 1 ≤ i ≤ n,
 
j6=i
es decir,
n n
! !
[ \ [
σ(A) ⊂ Fi Gi .
i=1 i=1
Definición 1.8 Sea f : C −→ C y A una matriz, diremos que f pertenece al espacio F(A)
si existe un entorno V de σ(A) sobre el cual f es analı́tica.
TEOREMA 1.6 (Teorema de la Aplicación Espectral) [11].

Si f ∈ F(A), entonces
f (σ(A)) = σ (f (A)) ,
donde f (σ(A)) = {f (λ) : λ ∈ σ(A)} .
LEMA 1.3 (Lema de Banach) [13]. Si P y Q son matrices de Cn×n y P es invertible,

−1
entonces si kP − Qk < P −1 se verifica que

Q es invertible , y P −1 − Q−1 < kP kkQkkP − Qk .

Observación 1.6 (Complemento de Schur) [4, pág. 92]. Sea A una matriz cuadrada
de orden n, con una partición en bloques

E F
A= ,
G H
donde la submatriz E es invertible y H cuadrada ( E y H pueden ser de tamaños distintos),

se llama complemento de Schur de E en A a la matriz
W = H − GE −1 F .
Ya que E es invertible, se tiene que

I E −1 F

I 0 E 0
A= , (1.7)
GE −1 I 0 H − GE −1 F 0 I
y puesto que el primer y el tercer factor del segundo miembro de (1.7) son invertibles se
deduce que A es invertible si y sólo si el complemento de Schur H − GE −1 F es invertible.
Observación 1.7 Dada una ecuación diferencial de 2o orden
Lu = a0 (x)u00 + a1 (x)u0 + a2 (x)u = f , 0 ≤ x ≤ 1,
se requieren las condiciones generales de la forma
α0,1 u(0) + α1,1 u0 (0) + β0,1 u(1) + β1,1 u0 (1) = γ1

, (1.8)
α0,2 u(0) + α1,2 u0 (0) + β0,2 u(1) + β1,2 u0 (1) = γ2
tales que

α0,1 α1,1 β0,1 β1,1
rango = 2, (1.9)
α0,2 α1,2 β0,2 β1,2
i.e., que ambas condiciones de contorno sean linealmente independientes.
Ahora trasladamos ésto a un sistema de ecuaciones en derivadas parciales de 2o orden con

coeficientes matriciales en las condiciones de contorno:

A1 u(0, t) + B1 ux (0, t) = 0
, (1.10)
A2 u(1, t) + B2 ux (1, t) = 0
donde Ak , Bk , k = 1, 2 son matrices cuadradas en Cm×m . Las condiciones de (1.10) son
un caso particular de las dadas en (1.8). Por tanto, en vista de la condición (1.9) se tiene
que verificar:

A1 B1
A = rango = 2m , (1.11)
A2 B2
i.e., A es una matriz invertible (ya que A es una matriz cuadrada de tamaño 2m × 2m) .
Los bloques de condiciones de contorno que aparecen en los problemas (1.1), (1.2) y (1.5)
son:

I 0 A1 B1 I 0
I) , II) , III) ,
C D A2 B2 B C
que por (1.11) tienen que ser invertibles. Por tanto, bajo la hipótesis de que A1 es invertible
y utilizando el complemento de Schur, véase la observación 1.6, obtenemos las siguientes
condiciones equivalentes:
I) invertible ↔ D invertible
II) invertible ↔ B2 − A2 A−1
1 B1 invertible
III) invertible ↔ C invertible .
10 1.2 Preliminares
Definición 1.9 Un polinomio p(λ), no nulo, se dice que es un polinomio anulador de A

si p(A) = 0.
Definición 1.10 Sea q(λ) el polinomio anulador de A de grado mı́nimo. Entonces a q(λ)
se le llama polinomio minimal de A.
TEOREMA 1.7 (Teorema de Cayley-Hamilton) [33, pág. 206].

Sea A una matriz cuadrada con polinomio caracterı́stico ϕ(λ), entonces ϕ(A) = 0. En otras
palabras, toda matriz cuadrada satisface su propia ecuación caracterı́stica.
Como consecuencia se tiene que el polinomio minimal de una matriz divide a su polinomio
caracterı́stico.
Para calcular el polinomio minimal de una matriz será suficiente tener en cuenta la defini-
ción 1.10 y el teorema 1.7.
Observación 1.8 El teorema 1.7 nos proporciona un método fácil para expresar cualquier
polinomio en A como otro, también en A, de grado menor que el del polinomio caracterı́stico
de A.
En efecto, si p(λ) es cualquier polinomio de grado mayor que n, entonces por el algoritmo
de la división, existen polinomios c(λ) y r(λ) no nulos tales que:
p(λ) = c(λ)ϕ(λ) + r(λ) (1.12)
donde ϕ(λ) es el polinomio caracterı́stico de A de grado n y r(λ) tiene grado menor que n.
Por tanto, aplicando el teorema 1.7 y utilizando (1.12) se tiene:
p(A) = c(A)ϕ(A) + r(A) = r(A) ,
i.e., el polinomio en A de grado mayor que n, p(A), se expresa como otro polinomio en A
de grado menor que n.
Notar que el mismo argumento también es válido si en lugar de considerar el polinomio
caracterı́stico de A se toma su polinomio minimal.
Definición 1.11 Se dice que los polinomios R1 (x), R2 (x), . . ., Rk (x) son primos 2 a 2
(coprimos), si Ri (x) y Rj (x) no tienen ningún divisor común, i 6= j. Se dice que son
primos entre sı́, si no tienen ningún divisor común a todos ellos.
Observación 1.9 En general, ser primos 2 a 2 implica ser primos entre sı́, pero ser primos
entre sı́ no implica ser primos 2 a 2.
TEOREMA 1.8 (Teorema de Bezout) [12]. Si Q1 (x) y Q2 (x) son dos polinomios pri-
mos entre sı́, entonces existen dos polinomios H1 (x) y H2 (x) de forma que se satisface la
relación:
1 = H1 (x)Q1 (x) + H2 (x)Q2 (x) .
TEOREMA 1.9 (Teorema de descomposición) [12, pág. 538].

Sea R(x) = R1 (x)R2 (x) · · · Rk (x) ∈ C[x] un polinomio factorizado, donde los factores
R1 (x), R2 (x), · · · , Rk (x) son primos 2 a 2. Consideremos el endomorfismo f ∈ L (Cn , Cn ).
Entonces:
Nuc R(f ) = Nuc R1 (f ) ⊕ Nuc R2 (f ) ⊕ · · · ⊕ Nuc Rk (f ) .
TEOREMA 1.10 [11, pág. 567]. Si P ∈ Cs×s , f (w) es una función holomorfa definida en
un conjunto abierto Ω del plano complejo y σ(P ) está en Ω, el cálculo funcional matricial
holomórfico define a f (P ) como una matriz que puede ser calculada como un polinomio en
P de grado más pequeño que el polinomio minimal de P .
COROLARIO 1.2 [40,pág. 76]. En particular, si P es invertible, entonces σ(P ) está en

iα

Dα = C ∼ Hα , Hα = −re : r ≥ 0 y considerando f (w) = logα (w) una rama del
√
logaritmo complejo asociado, holomórfico en Dα , entonces para w ∈ Dα , la función w =
exp 12 logα (w) es holomorfa y P = exp 12 logα (P ) es una raı́z cuadrada de P .

√
Observación 1.10 Si JP es la forma canónica
√ de Jordan de P , P = SJP S −1 y JP es
la raı́z cuadrada de JP , entonces Q = S JP S −1 es una raı́z cuadrada de P .
ESQUEMAS EN DIFERENCIAS FINITAS [45, pág. 13].
El método de diferencias finitas es uno de los métodos más utilizados para resolver proble-
mas de valores frontera. Dicho método consiste en reemplazar cada una de las derivadas de
la ecuación diferencial por el cociente incremental adecuado. A la ecuación discreta resul-
tante se le llama esquema en diferencias finitas.
Para poder realizar esta aproximación necesitamos definir una red (mallado) de puntos
en el plano (x, t). Tomando h (paso espacial) y k (paso temporal) números positivos, el
mallado estará formado por los puntos (xm , tn ) = (mh, nk) siendo m y n números enteros
arbitrarios.
La propiedad más básica que un esquema en diferencias finitas debe tener para ser útil es
la de ser un esquema convergente.
Definición 1.12 (Esquema en diferencias convergente) Sea P u = f una ecuación

en derivadas parciales siendo P un operador diferencial, y sea Pk,h v = f su esquema en
diferencias finitas. Se dice que el esquema Pk,h v = f es convergente si su solución tiende a
la solución de la ecuación en derivadas parciales cuando h, k −→ 0.
Pero como probar la convergencia de un esquema en general no es fácil, se utilizan dos

conceptos relacionados con la convergencia que son más fáciles de analizar: la consistencia
y la estabilidad. Utilizaremos que la consistencia y la estabilidad implican la convergencia
(teorema de Lax ).
Definición 1.13 (Esquema en diferencias consistente)

Dada una ecuación en derivadas parciales P u = f y su esquema en diferencias finitas
Pk,h v = f , diremos que el esquema en diferencias finitas es consistente respecto a la ecua-
ción en derivadas parciales si para cualquier función φ(x, t) suave (suficientemente deriva-
ble para el contexto) se satisface
P φ − Pk,h φ −→ 0 , cuando k, h → 0,
donde la convergencia es puntual en cada punto del mallado.

En otras palabras, un esquema Pk,h v = f es consistente respecto a la ecuación P u = f si
la solución de la ecuación, si ésta es suave, es una solución aproximada del esquema en
diferencias.
12 1.2 Preliminares
Definición 1.14 (Esquema en diferencias estable) Sea v(mh, nk) una solución del es-
quema Pk,h v = f . Diremos que el esquema es estable si para un h0 (paso espacial) fijado,
la solución v(mh0 , nk) permanece acotada cuando n → ∞.
Definición 1.15 (Estabilidad en el sentido de estación fija) [43, pág. 109] Conside-
remos una solución v(mh, nk) del esquema Pk,h v = f . Diremos que el esquema es estable
m
si dado un punto del mallado (X, T ), donde X = mh0 = M siendo h0 fijo y T = Jk finito,
entonces la solución v(mh0 , Jk) permanece acotada a medida que el paso temporal k decrece
(k → 0), i.e., a medida que n crece, dependiendo del valor de k, hasta tomar el valor J
donde se alcanza el tiempo T .
El estudio de la estabilidad numérica es relevante, dado que las aproximaciones numéricas

siempre están sometidas a pequeñas perturbaciones: los errores de redondeo. Estos errores
son debidos a que en los cálculos prolongados es probable que se realicen muchos redondeos,
donde cada uno de ellos desempeña el papel de un error de entrada para el resto del cálculo
y cada uno tiene un efecto sobre la consiguiente salida. La obtención de soluciones estables
nos garantizará el tener controlados estos errores.
En esta memoria el estudio de la estabilidad se ha realizado diferenciando entre una esta-

bilidad uniforme en el tiempo, es decir, independientemente de t (véase la definición 1.14)
y una estabilidad puntual, es decir, relativa a un punto concreto del mallado (véase la
definición 1.15). En ambas definiciones se ha considerado que el paso espacial h está fijo.
A continuación comentaremos varios teoremas y definiciones de los problemas de Sturm-

Liouville discretos, que serán de gran utilidad en este trabajo. Todos los resultados y sus
demostraciones pueden encontrarse en [1, cap. 11].
PROBLEMAS DE STURM-LIOUVILLE DISCRETOS
Obviamente, los problemas de valores frontera lineales homogéneos pueden tener soluciones
no triviales. Si los coeficientes de la ecuación en diferencia y/o los de las condiciones de
frontera dependen de un parámetro, entonces uno de los problemas pioneros de la fı́sica
matemática es determinar el/los valor(es) del parámetro para el cual tales soluciones no tri-
viales existan. Esos valores especiales del parámetro se llaman autovalores o valores propios
y las correspondientes soluciones no triviales se llaman autofunciones o funciones propias.
El problema de valores frontera formado por la ecuación en diferencias
∆ (p(k − 1)∆u(k − 1)) + q(k)u(k) + λr(k)u(k) = 0 , k ∈ N (1, K) , (1.13)
y las condiciones de frontera
u(0) = αu(1) , u(K + 1) = βu(K) , (1.14)
se llama problema de Sturm-Liouville discreto (P.S.L.D.). Donde el operador de

diferencia, ∆, viene dado por
∆p(k − 1) = p(k) − p(k − 1) , y ∆2 p(k) = ∆ (∆p(k)) .

En la ecuación de diferencia (1.13), λ es el parámetro, y las funciones p, q y r están definidas

sobre los conjuntos N (0, K), N (1, K) y N (1, K) respectivamente, y p(k) > 0, k ∈ N (0, K),
r(k) > 0, k ∈ N (1, K). Los conjuntos N (0, K) y N (1, K) están definidos como sigue
N (0, K) = {a ∈ N / 0 ≤ a ≤ K} ,
N (1, K) = {b ∈ N / 1 ≤ b ≤ K} .
En las condiciones de frontera (1.14), α y β son constantes conocidas.
Los siguientes resultados, en los cuales se asume tácitamente la existencia de los valores
propios del P.S.L.D. (1.13)-(1.14), son fundamentales.
TEOREMA 1.11 Los autovalores del P.S.L.D. (1.13)-(1.14) son simples, i.e., si λ es
un autovalor de (1.13)-(1.14) y φ1 (k) y φ2 (k) son sus respectivas autofunciones, entonces
φ1 (k) y φ2 (k) son linealmente dependientes en N (0, K + 1).
Definición 1.16 El conjunto de funciones {φm (k) / m = 1, 2, . . .} cada una de las cuales
está definida en N = N ∪ {0} se dice ortogonal en N con respecto a la función no negativa
r(k), k ∈ N si X
r(l)φµ (l)φν (l) = 0 , ∀µ 6= ν .
l∈N
La función r(k) se llama función peso.
TEOREMA 1.12 Sean λm ; m = 1, 2, . . ., los autovalores del P.S.L.D. (1.13)-(1.14), y

φm (k); m = 1, 2, . . ., las correspondientes autofunciones. Entonces, el conjunto {φm (k) / m = 1, 2, . . .}
es ortogonal en N (1, K) con respecto a la función peso r(k).
TEOREMA 1.13 Sean λ1 y λ2 dos autovalores del P.S.L.D. (1.13)-(1.14) y sean φ1 (k)
y φ2 (k) las correspondientes autofunciones. Entonces φ1 (k) y φ2 (k) son linealmente depen-
dientes en N (0, K + 1) si y sólo si λ1 = λ2 .
TEOREMA 1.14 Para el P.S.L.D. (1.13)-(1.14), los autovalores son reales.
FORMULACIÓN MATRICIAL DE PROBLEMAS DE STURM-LIOUVILLE

DISCRETOS
Al desarrollar la ecuación en diferencia (1.13), utilizando la definición del operador de

diferencia ∆, obtenemos para k ∈ N (1, K):
p(k)u(k + 1) − (p(k) + p(k − 1)) u(k) + (q(k) + λr(k)) u(k)

(1.15)
+p(k − 1)u(k − 1) = 0 .
Si denotamos s(k) = p(k) + p(k − 1) − q(k), con k ∈ N (1, K), entonces la ecuación (1.15)
se puede reescribir ası́
− p(k − 1)u(k − 1) + s(k)u(k) − p(k)u(k + 1) = λr(k)u(k) . (1.16)

k ∈ N (1, K)
14 1.2 Preliminares
Por tanto, para k = 1, . . . , K, tenemos las siguientes ecuaciones:
− p(0)u(0) + s(1)u(1) − p(1)u(2) = λr(1)u(1) , (1.17)
y
− p(K − 1)u(K − 1) + s(K)u(K) − p(K)u(K + 1) = λr(K)u(K) , (1.18)
las cuales en vista de las condiciones de contorno (1.14) toman la forma
s(1)u(1) − p(1)u(2) = λr(1)u(1) , (1.19)
− p(K − 1)u(K − 1) + s(K)u(K) = λr(K)u(K) , (1.20)

donde s(1) = s(1) − αp(0) y s(K) = s(K) − βp(K).
Las K ecuaciones: (1.19), (1.16) para k ∈ N (2, K − 1) y (1.20) pueden ser escritas como
un problema de autovalores matricial
A
eu = λ R u , (1.21)
e es la matriz real, K × K, simétrica y tridiagonal de la forma HK (F, G), con
donde A
F = (s(1), s(2), . . . , s(K − 1), s(K)) ,
y
G = (−p(1), −p(2), . . . , −p(K − 1)) ;
R es la matriz diagonal K × K definida como R = diag(r(1), . . . , r(K)), y la incógnita
u = (u(1), u(2), . . . , u(K)).
Ya que p(k) > 0, k ∈ N (0, K) y r(k) > 0, k ∈ N (1, K) se sigue que:
(i) el P.S.L.D (1.13)-(1.14) tiene exactamente K autovalores reales λm , 1 ≤ m ≤ K, los

cuales son distintos;
(ii) correspondiendo a cada autovalor, λm , existe un autofunción φm (k), k ∈ N (1, K).

Esas autofunciones φm (k), 1 ≤ m ≤ K son mútuamente ortogonales con respecto a la
función peso r(k). En particular, esas autofunciones son linealmente independientes
en N (1, K).
SERIES DE FOURIER DISCRETAS
Sean a, b ∈ Z y {φm (k) / a ≤ m ≤ b} un conjunto ortogonal de funciones en N (a, b) con

respecto a la función peso positiva r(k), k ∈ N (a, b). Ya que la ortogonalidad de esas
funciones φm (k), a ≤ m ≤ b, en particular, implica su independencia lineal en N (a, b),
entonces cualquier función u(k), k ∈ N (a, b) puede ser expresada como una combinación
lineal de φm (k), a ≤ m ≤ b, i.e.,
b
X
u(k) = cm φm (k) , k ∈ N (a, b) , (1.22)
m=a
donde las constantes cm , a ≤ m ≤ b, pueden ser determinadas como sigue:

(1) multiplicando ambos lados de (1.22) por r(k)φn (k), a ≤ n ≤ b;
(2) sumando el resultado desde k = a hasta k = b;
(3) utilizando la ortogonalidad de las funciones φm (k), a ≤ m ≤ b, en N (a, b).
Por tanto, obtenemos:

b b b
!
X X X
r(k)φn (k)u(k) = cm r(k)φn (k)φm (k)
k=a m=a k=a
b
X
= cn r(k)φ2n (k) ,
k=a
y entonces
b
P
r(k)φm (k)u(k)
k=a
cm = b
, a ≤ m ≤ b. (1.23)
r(k)φ2m (k)
P
k=a
En particular, si las funciones φm (k), a ≤ m ≤ b son ortonormales, i.e., para cada m,

b
r(k)φ2m (k) = 1, entonces las constantes cm se simplifican a
P
k=a
b
X
cm = r(k)φm (k)u(k) , a ≤ m ≤ b.
k=a
Definición 1.17 La relación (1.22) se llama la serie de Fourier discreta, y las constantes
cm definidas en (1.23) son los correspondientes coeficientes de Fourier discretos.
Capı́tulo 2
Sistemas parabólicos
2.1. Caso particular

0.4pt0.4pt 0pt0.4pt
2.1.1. Introducción
Este capı́tulo trata sobre la construcción de soluciones numéricas convergentes de problemas
mixtos con condiciones de contorno acopladas del tipo:
ut (x, t) − Auxx (x, t) − Bu(x, t) = 0, 0 < x < 1, t > 0, (2.1)

u(0, t) = 0, t > 0, (2.2)
Cu(1, t) + Dux (1, t) = 0, t > 0, (2.3)
u(x, 0) = F (x), 0 ≤ x ≤ 1, (2.4)
donde u = (u1 , u2 , . . . , us )T y F (x) = (f1 , . . . , fs )T son vectores s-dimensionales, elementos

de Cs y A, B, C y D son s × s matrices complejas, elementos de Cs×s verificando ciertas
propiedades a determinar. La solución numérica propuesta es la solución exacta del esque-
ma en diferencias obtenido al discretizar el problema continuo (2.1)–(2.4). Dicha solución
se obtiene utilizando un método discreto de separación de variables que evita el cálculo de
los sistemas algebraicos de gran tamaño con bloques matriciales que aparecen en el método
de diferencias estándar. Para el caso D = 0 y C invertible, el problema (2.1)–(2.4) ha sido
tratado en [37] utilizando un método de separación de variables matricial. El acoplamiento
en la condición de contorno (2.3) y la ley de conmutatividad entre matrices hace difı́cil
el tratamiento analı́tico del problema (2.1)–(2.4) utilizando un método de separación de
variables.
En el apartado 2.1.2, estudiaremos la discretización del problema mixto (2.1)–(2.4) utili-

zando aproximaciones en diferencias centradas para la derivada de segundo orden uxx , y
aproximaciones en diferencias progresivas para ut . En el apartado 2.1.3, trataremos la exis-
tencia y construcción de soluciones no triviales para el problema de contorno (2.1)–(2.3)
utilizando un método de separación de variables discreto. En el apartado 2.1.4 se compro-
bará que el esquema en diferencias finitas es consistente respecto a la ecuación (2.1) y se
17
18 2.1 Caso particular
buscarán condiciones suficientes para garantizar que las soluciones construidas sean esta-
bles. De esta forma aseguraremos la convergencia de las soluciones. En el apartado 2.1.5,
se construirán soluciones para el problema mixto (2.1)–(2.4) discretizado. En el apartado
2.1.6, utilizaremos un método de proyecciones que nos permitirá extender los resultados
del apartado 2.1.5 a clases más amplias de funciones de valores iniciales F (x). Finalmente,
en el apartado 2.1.7 desarrollaremos algunos ejemplos para aplicar los resultados teóricos
obtenidos.
2.1.2. La discretización del problema

Para la discretización del dominio [0, 1] × [0, +∞[, se considera una discretización de [0, 1]
1
tomando un número entero M > 0 y un paso espacial, h = M que consideraremos fijo para
simplificar, i.e., se construye el conjunto de puntos:
{xm = mh, m = 0, 1, . . . , M, } ,
por otra parte, se toma un paso temporal k > 0 (también fijo por la misma razón), y se
discretiza el tiempo, [0, +∞[, considerando la sucesión de puntos:
{tn = nk, n = 0, 1, . . .} .
Con esto hemos construido una red rectangular del dominio [0, 1]×[0, +∞[, formada por los
puntos (xm , tn ) = (mh, nk). En cada uno de ellos utilizaremos un esquema en diferencias
finitas progresivas en el que se reemplazarán las derivadas parciales ut y uxx , de la ecuación
(2.1), convirtiendo el problema continuo en uno discreto. Para reemplazar ut desarrollamos
u(mh, (n + 1)k) en serie de Taylor, en t, alrededor del punto (mh, nk) obteniendo:
u (mh, (n + 1)k) = u(mh, nk) + kut (mh, nk) + O(k 2 ) ,

despejando ut (mh, nk):
u (mh, (n + 1)k) − u(mh, nk)

ut (mh, nk) = + O (k) . (2.5)
k
De la misma forma que para ut , podemos reemplazar uxx por una aproximación en diferen-
cias finitas centradas desarrollando u ((m + 1)h, nk) y u ((m − 1)h, nk) en serie de Taylor,
en x, alrededor del punto (mh, nk):
u ((m + 1)h, nk) =

h2 h3
u(mh, nk) + hux (mh, nk) + uxx (mh, nk) + uxxx (mh, nk) + O(h4 ) ,
2! 3!
u ((m − 1)h, nk) =

h2 h3
u(mh, nk) − hux (mh, nk) + uxx (mh, nk) − uxxx (mh, nk) + O(h4 ) .
2! 3!
Sumando estas dos últimas expresiones y despejando el término uxx tenemos:
Capı́tulo 2. Sistemas parabólicos 19
u ((m + 1)h, nk) − 2u(mh, nk) + u ((m − 1)h, nk)

uxx (mh, nk) = + O(h2 ) . (2.6)
h2
Entonces sustituyendo en el problema (2.1)–(2.4) las expresiones (2.5) y (2.6), y denotando

U (m, n) = u(mh, nk) obtenemos el siguiente problema mixto discreto
rA [U (m + 1, n) + U (m − 1, n)]

rB
+ I + 2 − 2rA U (m, n) − U (m, n + 1) = 0 , 0 < m < M, n > 0 (2.7)
M
U (0, n) = 0 , n > 0, (2.8)

CU (M, n) + M D [U (M, n) − U (M − 1, n)] = 0 , n > 0, (2.9)
U (m, 0) = f (m) , 0 ≤ m ≤ M, (2.10)
donde
k 1
r= , h= , f (m) = F (mh) . (2.11)
h2 M
2.1.3. El problema de contorno discreto

En este apartado buscamos soluciones no triviales del problema de contorno (2.7)–(2.9) de
la forma
U (m, n) = G(n)H(m), G(n) ∈ Cs×s , H(m) ∈ Cs . (2.12)

Sustituyendo (2.12) en (2.7) obtenemos
rAG(n) [H(m + 1) + H(m − 1)] +

rB
I + 2 − 2rA G(n)H(m) − G(n + 1)H(m) = 0. (2.13)
M
Tomando un número real ρ, la ecuación (2.13) puede escribirse de la forma

2r + ρ
rAG(n) H(m + 1) − H(m) + H(m − 1)
r
(2.14)
rB
+ I + 2 + ρA G(n) − G(n + 1) H(m) = 0 ,
M
y (2.14) se cumple si {H(m)}, {G(n)} satisfacen

rB
G(n + 1) − I + 2 + ρA G(n) = 0 , n > 0 , (2.15)
M

2r + ρ
H(m + 1) − H(m) + H(m − 1) = 0 , 1 ≤ m ≤ M − 1 . (2.16)
r
Observar que si ρ satisface
− 4r < ρ < 0 , (2.17)
2
2r+ρ
entonces 2r ≤ 1. Por tanto la ecuación algebraica

2 2r + ρ
z − z+1=0, (2.18)
r
tiene dos soluciones diferentes dadas por

r 2 
z0 = 2r + i 1 − 2r+ρ
2r+ρ
2r = eiθ , 


r 
2 , (2.19)
2r+ρ 2r+ρ
z1 = 2r −i 1− 2r = e−iθ , 



2r+ρ θ
, ρ = −4r sen2 , i2 = −1

donde 0 < θ < π , cos θ = 2r 2
y observamos que
z0n = cos(nθ) + i sen(nθ); z1n = cos(nθ) − i sen(nθ) .
Ya que la ecuación vectorial (2.16) tiene coeficientes escalares, sus soluciones pueden escri-
birse de la forma
H(m) = cos(mθ) c + sen(mθ) d , c, d ∈ Cs , 1 ≤ m ≤ M − 1. (2.20)
Por (2.8) y (2.12), teniendo en cuenta que {G(n)} no puede ser idénticamente cero, se sigue
que H(0) = 0 y por (2.20) obtenemos
H(m) = sen(mθ)d , d ∈ Cs , 1 ≤ m ≤ M − 1 . (2.21)

La solución de (2.15) verificando G(0) = I viene dada por
n
rB
G(n) = I + 2 + ρA , n > 0. (2.22)
M
Sustituyendo (2.21) y (2.22) en (2.12) obtenemos

n
rB
U (m, n) = I + 2 + ρA sen(mθ) d , d ∈ Cs , n > 0 , 1 ≤ m ≤ M − 1 ,
M
e imponiendo la condición de contorno (2.9), resulta que el vector d debe satisfacer
{C sen(M θ) + M D [sen(M θ) − sen ((M − 1)θ)]} (I + ρ(A + wB))n d = 0 , (2.23)

donde
r
w= 6= 0 , (2.24)
ρM 2
y ρ viene dado por (2.19). Ya que por el teorema de Cayley-Hamilton, véase el teorema 1.7,
las potencias (A + wB)n pueden expresarse en términos de I, A + wB, (A + wB)2 , . . .,
(A + wB)p−1 , donde p es el grado del polinomio minimal de la matriz A + wB, se tiene que
la condición (2.23) es equivalente a la condición
{C sen(M θ) + M D [sen(M θ) − sen ((M − 1)θ)]} (A + wB)n d = 0 , (2.25)

0 ≤ n < p.
Como estamos interesados en soluciones no nulas, el vector d debe ser no nulo y por tanto
es suficiente que
L(θ) = C sen(M θ) + M D [sen(M θ) − sen ((M − 1)θ)] es singular , (2.26)

0 < θ < π.
`π
Observar que si θ` = M para 1 ≤ ` ≤ M − 1, entonces sen(M θ) = 0 y se obtiene que en
`π
este caso la matriz L(θ` ) = M D(−1)` sen M

es invertible si D es invertible, puesto que
`π
sen( M ) 6= 0. Por tanto bajo la hipótesis
D es invertible , (2.27)
véase la observación 1.7, los valores de θ para los cuales se verifica (2.26) deben satisfacer
que el sen(M θ) 6= 0 y podemos escribir (2.26) de la forma
M [sen(M θ) − sen ((M − 1)θ)]

D−1 C + I es singular. (2.28)
sen(M θ)
Asumamos que
Existe µ ∈ σ(−D−1 C) ∩ R . (2.29)
La condición (2.28) se cumple si existen θ para la ecuación escalar
M [sen(M θ) − sen ((M − 1)θ)]
= µ. (2.30)
sen(M θ)
La ecuación (2.30) es equivalente a
µ
cos θ − 1 − M
cot(M θ) = . (2.31)
sen θ
i h
(`−1)π `π
Observamos que para cada ` con 1 ≤ ` ≤ M −1, en el intervalo J` = M ,M se satisface
lı́m cot(M θ) = +∞ ; lı́m cot(M θ) = −∞ ; cot(M θ) decrece en J` , (2.32)

(`−1) + ` −
θ→ π θ→ M π
M
d
porque dθ (cot(M θ)) = − sen2M(M θ)
< 0. Para buscar las raı́ces de la ecuación (2.31) distin-
guiremos los siguientes casos sobre el valor propio real µ de la matriz −D−1 C.
CASO 1. µ<0
Para este caso la función eM (θ) definida por

µ
cos θ − 1 − M
eM (θ) = , (2.33)
sen θ

M
es continua, negativa y tiene un máximo en θ = arc cos M −µ a lo largo de ]0, π[.
CASO 2. µ=0
Si C es singular, µ toma el valor 0 y (2.30) es equivalente a la ecuación
sen(M θ) = sen((M − 1)θ) ,
que tiene M − 1 soluciones dadas por θ` = (2`−1)π

2M −1 para 1 ≤ ` ≤ M − 1. Observamos
M −1
entonces que la matriz
L(θ singular para los valores de {θ` }`=1 obtenidos, puesto que
` ) es
`−M
C es singular, sen 2M −1 π 6= 0 para 1 ≤ ` ≤ M − 1 y L(θ` ) toma la forma

(2` − 1) (2` − 1)
L(θ` ) = 2M D sen π cos π
2(2M − 1) 2

2` − 1
+ C sen M π
2M − 1

`−M
= (−1)` C sen π , 1 ≤ ` ≤ M − 1.
2M − 1
Por tanto, por lo visto en los Casos 1 y 2, se obtiene que si µ ≤ 0 existe una única raı́z,
θ` ∈ J` de la ecuación (2.31) para 1 ≤ ` ≤ M − 1, es decir, existen M − 1 soluciones θ` de
la ecuación (2.31) en ]0, π[.
CASO 3. µ ∈ ]0, +∞[

π
En el intervalo J1 = 0, M no podemos garantizar analı́ticamente laiexistencia de
h solución
para la ecuación (2.31). En cambio para el resto de intervalos, J` = (`−1) `
M π, M π , 2 ≤ ` ≤
M − 1, existe una única solución θ` ∈ J` de la ecuación (2.31) porque la función eM (θ) es
continua, decreciente y acotada en cada J` ( ya que lı́mθ→ (`−1) π+ eM (θ) y lı́mθ→ ` π− eM (θ)
M M
son números reales).
Por tanto, si µ ∈ ]0, +∞[ existe una única raı́z, θ` ∈ J` de la ecuación (2.31) para 2 ≤ ` ≤
M − 1, es decir, existen M − 2 soluciones θ` de la ecuación (2.31) en ]0, π[.
Sea µ un valor propio real de −D−1 C y θ` la solución de (2.31) en cada J` . Recordemos

que si µ ≤ 0 entonces 1 ≤ ` ≤ M − 1 y que para µ ∈ ]0, +∞[ se tiene 2 ≤ ` ≤ M − 1.
e (µ, θ` ) ∈ Csp` ×s definidas por

Introducimos las matrices G(µ) ∈ Cs×s y G
D−1 L(θ` )
G(µ) = = D−1 C + µI , (2.34)
sen(M θ` )
 
G(µ)
 G(µ) (A + w` B) 

G(µ) (A + w` B)2
 −1
G
e (µ, θ` ) =  , w` = , (2.35)


 ..  4M 2 sen2 θ`
 .  2
G(µ) (A + w` B)p` −1
donde p` es el grado del polinomio minimal de la matriz A + w` B. Se observa que utilizando
la matriz G
e (µ, θ` ), la condición (2.25) toma la forma

θ` 2M − 1
C sen (M θ` ) + 2M D sen cos θ` (A + w` B)n d` = 0,
2 2
0 ≤ n < p` ,
es decir,
L(θ` ) (A + w` B)n d` = 0 , 0 ≤ n < p` , (2.36)
que es equivalente a la condición
G
e (µ, θ` ) d` = 0, d ` ∈ Cs . (2.37)
La ecuación (2.37) tiene soluciones no triviales si y sólo si

rango G e (µ, θ` ) < s , (2.38)
y por el teorema de Mitra, véase el teorema 1.1, bajo la hipótesis (2.38) el conjunto de
soluciones de la ecuación algebraica (2.37) viene dado por

d` = I − G e (µ, θ` )† G
e (µ, θ` ) S` , S` ∈ Cs ∼ {0}. (2.39)
Resumiendo, bajo las hipótesis (2.27) y (2.29), si θ` son las soluciones de (2.31) y G
e (µ, θ` ) la
matriz definida por (2.35), el conjunto de soluciones no triviales del problema de contorno
(2.7)–(2.9) viene dado por
n 
U` (m, n) = I − r 4 sen2 θ2` A − B
M2
sen(mθ` ) d` 



1 ≤ m ≤ M − 1, n > 0

, (2.40)


e (µ, θ` )† G

d` = I − G e (µ, θ` ) S` , S` ∈ Cs ∼ {0}  
donde 1 ≤ ` ≤ M − 1 si µ ≤ 0 y 2 ≤ ` ≤ M − 1 si µ ∈ ]0, +∞[.
Supongamos ahora que ρ = 0, lo cual corresponde a θ = 0 en (2.19). En este caso la

ecuación (2.18) toma la forma:
z 2 − 2z + 1 = (z − 1)2 = 0,
y el conjunto de soluciones de la ecuación (2.16) viene dado por:
H(m) = 1m c + m 1m d1 = c + m d1 , c, d1 ∈ Cs , 1 ≤ m ≤ M − 1. (2.41)
Sabı́amos que era necesario que H(0) = 0 para no tener soluciones triviales del problema
de frontera discreto (2.7)–(2.9). Por tanto, de (2.16) para m = 1 y ρ = 0, y de (2.41) para
m = 1, 2, se sigue que c = 0. Por otro lado, de (2.16) para m = M − 1 y ρ = 0, y de (2.41)
para m = M − 2, M − 1, obtenemos H(M ) = M d1 . De esta forma se llega a que
H(m) = m d1 , d 1 ∈ Cs , 1 ≤ m ≤ M. (2.42)
Como ρ = 0, de (2.15) se sigue que:

r n
G(n) = I + 2 B , n > 0. (2.43)
M
Si sustituimos (2.42) y (2.43) en la condición frontera (2.9) llegamos a la siguiente expresión:
r n
C I + 2 B M d1
M h r n r n i
+ MD I+ B M d 1 − I + B (M − 1) d 1 = 0, n > 0,
M2 M2
o equivalentemente
r n
(C + D) I + 2 B d1 = 0, d1 ∈ Cs , n > 0. (2.44)
M
Por el teorema de Cayley-Hamilton, véase el teorema 1.7, si t es el grado del polinomio
minimal de B entonces para n ≥ t, las potencias B n pueden ser expresadas en términos de
las matrices I, B, B 2 , . . . , B t−1 . Como Mr 2 = k 6= 0, para que se cumpla (2.44) es suficiente
considerar:
(C + D)B n d1 = 0, d1 ∈ Cs , 0 ≤ n < t. (2.45)
Como D es invertible, la ecuación (2.45) equivale a
D−1 C + I B n d1 = 0, d1 ∈ Cs ,

0 ≤ n < t. (2.46)
Puesto que d1 ∈ Cs debe ser no nulo, para que se verifique (2.46) es suficiente que se
cumpla:
D−1 C + I singular . (2.47)
La condición (2.47) es equivalente a que µ = 1 ∈ σ −D−1 C . Supongamos que ∃ µ = 1 ∈

σ(−D−1 C). A fin de determinar el vector d1 , introducimos la matriz R(1)

e ∈ Cst×s definida
por:
 
G(1)
 G(1)B 
 
2
R(1) =  G(1)B , (2.48)
e  
 .. 
 . 
G(1)B t−1
donde la matriz G(1) ha sido definida en (2.34). Observamos que utilizando R(1),
e la con-
dición (2.46) es equivalente a :
R(1)
e d1 = 0, d 1 ∈ Cs . (2.49)
La ecuación (2.49) tiene una solución no nula

d1 ∈ Cs ⇐⇒ rango R(1)
e < s,
y bajo esta condición, por el teorema de Mitra, véase el teorema 1.1, el conjunto solución
de (2.49) viene dado por:

e † R(1)
d1 = I − R(1) e S1 , S1 ∈ Cs ∼ {0}. (2.50)
Por tanto, para µ = 1 obtenemos la solución no trivial del problema de frontera (2.7)–(2.9):
n
U (m, n) = m I + Mr 2 B d1 , 1 ≤ m ≤ M − 1, n > 0 


. (2.51)
e † R(1)
d1 = I − R(1) e S1 , S1 ∈ Cs ∼ {0}


Resumiendo, hemos demostrado el siguiente resultado.
TEOREMA 2.1 Consideremos el problema de frontera (2.7)–(2.9) y asumamos la condi-

ción (2.27). Sea µ un número real verificando (2.29). Sea G(µ) la matriz en Cs×s definida
por (2.34). Si M > 0, r = hk2 , sean w` , Ge (µ, θ` ) y R(1)
e definidas por (2.35) y (2.48) res-
pectivamente, donde θ` es una solución de (2.31). Sean p` y t los grados de los polinomios
minimales de las matrices A + w` B y B respectivamente.
(i) Si µ ∈ R, existen valores θ ∈ ]0, π[ tales que la matriz L(θ) definida por (2.26) es
singular. Bajo la hipótesis (2.26)
el problema
(2.7)–(2.9) tiene soluciones no triviales
de la forma (2.12) si rango G e (µ, θ) < s.

2`−1
(ii) Si C es singular, se toma el valor propio µ = 0 entonces θ` = 2M −1 π para 1 ≤

` ≤ M − 1 son soluciones de (2.31). Si rango Ge (0, θ` ) < s entonces {U` (m, n)}M −1
`=1
dadas por (2.40) define M − 1 soluciones no triviales del problema de contorno (2.7)–
(2.9).
(iii) Si µ < 0, entonces

i para h1 ≤ ` ≤ M − 1 existe una solución θ` de la ecuación
(`−1)π `π
(2.31) en J` = M , M para la cual la matriz L(θ` ) es singular. Bajo la hipóte-

sis rango Ge (µ, θ` ) < s, la sucesión {U` (m, n)}M −1 dada por (2.40) define M − 1
`=1
soluciones no triviales del problema de contorno (2.7)–(2.9).
(iv) Si µ = 1, entonces i parah 2 ≤ ` ≤ M − 1 existe una solución θ` de la ecuación

(`−1)π `π
(2.31) en J` = , M para la cual la matriz L(θ` ) es singular. Bajo las hipótesis
M
rango G e (1, θ` ) < s y rango R(1) e < s, existen dos conjuntos de soluciones no
triviales del problema de contorno (2.7)–(2.9) dadas por (2.40) y (2.51).
(v) Si µ ∈ ]0, +∞[∼ {1}, entonces

i parah 2 ≤ ` ≤ M − 1 existe una solución θ` de la
(`−1)π `π
ecuación (2.31) en J` = , M para la cual la matriz L(θ` ) es singular. Bajo
M
la hipótesis rango Ge (µ, θ` ) < s, la sucesión {U` (m, n)}M −1 dada por (2.40) define
`=2
M − 2 soluciones no triviales del problema de contorno (2.7)–(2.9).
2.1.4. Consistencia y estabilidad de las soluciones

Veamos que el esquema en diferencias (2.7) es consistente respecto a la ecuación (2.1) en
el sentido de la definición 1.13. Denotamos
Λ[u] = ut (x, t) − Auxx (x, t) − Bu(x, t) , (2.52)

u (mh, (n + 1)k) − u(mh, nk)
Λk,h [u] = − Bu(mh, nk)
k (2.53)
[u ((m + 1)h, nk) − 2u(mh, nk) + u ((m − 1)h, nk)]
−A .
h2
Sea φ(x, t) una función suave, i.e. suficientemente derivable para el contexto, y sea Φ(m, n) =
φ(mh, nk). Considerando el desarrollo de Taylor de Φ para el punto (m, n + 1) con respecto
a la variable t, y para los puntos (m + 1, n), (m − 1, n) respecto de x, se sigue
2 3

Φ(m, n+1) = Φ(m, n)+kΦt(m, n)+ k2! Φtt(m, n)+ k3! Φttt(m, n)+ O(k 4 )


2 3
Φ(m+1, n) = Φ(m, n)+hΦx(m, n)+ h2! Φxx(m, n)+ h3! Φxxx(m, n)+ O(h4 ) (2.54)



2 3
Φ(m−1, n) = Φ(m, n)−hΦx(m, n)+ h2! Φxx(m, n)− h3! Φxxx(m, n)+ O(h4 )
Sustituyendo (2.54) en Λk,h [φ], dada por (2.53), se sigue que
k k2
Λk,h [φ] = Φt (m, n) + Φtt (m, n) + Φttt (m, n)
2 6
h2

−A Φxx (m, n) + Φxxxx (m, n)
12
−BΦ(m, n) + O(k 3 ) + O(h4 ) . (2.55)
Teniendo en cuenta la expresión de Λ[φ], dada por (2.52), y la expresión (2.55) se sigue que
k h2
Λ[φ] − Λk,h [φ] = − Φtt (m, n) + A Φxxxx (m, n) + O(k 2 + h4 ) . (2.56)
2 12
Tomando normas en (2.56) obtenemos
k Λ[φ] − Λk,h [φ] k ≤ O(k) + kAkO(h2 ) −→ 0 , cuando k, h → 0 . (2.57)

Por tanto el esquema en diferencias (2.7) es consistente respecto a la ecuación (2.1).
Estudiemos ahora la estabilidad, véase la definición 1.14, de las soluciones construidas en

el teorema 2.1, a fin de tener controlados los errores de redondeo que pueden aparecer
en el cómputo de dichas soluciones a medida que n crezca. A continuación buscaremos
condiciones suficientes para garantizar la estabilidad de las expresiones
n
2 θ` 1
I − r 4 sen A − 2B , n>0 (2.58)
2 M
r n
I + 2B , n > 0 (2.59)
M
cuando n → ∞, es decir para garantizar que las expresiones matriciales
r
−B + α2 A ,

I− 2
(2.60)
M
r
I+ B, (2.61)
M2

donde se ha tomado α = 2M sen θ2` , sean matrices convergentes, véase el lema 1.1.
Comenzaremos nuestro estudio con la expresión matricial dada en (2.60).
Por el teorema de la aplicación espectral, véase el teorema 1.6, se sabe que

r n r o
σ I − 2 −B + α2 A = 1 − 2 z : z ∈ σ −B + α2 A . (2.62)
M M
Si z está en σ −B + α2 A entonces

r 2 r 2 r 2 r2 2r
1 − 2 z = 1 − 2 Re(z) + − 2 Im(z) = 1 + 4 |z|2 − 2 Re(z) . (2.63)

M M M M M
Por (2.62), (2.63) la matriz dada en (2.60) es convergente si
r 2 Re(z)
para todo z ∈ σ −B + α2 A .

Re(z) > 0 y 2
< 2
(2.64)
M |z|
Puesto que A y B son matrices complejas en Cs×s se pueden descomponer de la siguiente

forma, véase la definición 1.7:
A = A1 + i A 2 , B = B1 + i B2 ,
donde Ak , Bk , para k = 1, 2, son las matrices hermı́ticas definidas por
B + BH B − BH A + AH A − AH
B1 = ; B2 = ; A1 = ; A2 = , (2.65)
2 2i 2 2i
entonces
−B + α2 A = −B1 + α2 A1 + i −B2 + α2 A2 .

(2.66)
Ya que α2 Ak − Bk , k = 1, 2, son matrices hermı́ticas podemos aplicar el teorema de Bro-

mwich, véase el teorema 1.3, obteniendo
λmı́n −B1 + α2 A1 ≤ Re(z) ≤ λmáx −B1 + α2 A1



2
2

λmı́n −B2 + α A2 ≤ Im(z) ≤ λmáx −B2 + α A2  . (2.67)

z ∈ σ −B + α2 A

Por el teorema 1.4 se sigue que
λmı́n −B1 + α2 A1 ≥ λmı́n (−B1 ) + α2 λmı́n (A1 )

= α2 λmı́n (A1 ) − λmáx (B1 ) , (2.68)
y por (2.67), (2.68) se obtiene
mı́n Re(z) : z ∈ σ −B + α2 A ≥ α2 λmı́n (A1 ) − λmáx (B1 ) .

(2.69)
Bajo las hipótesis
x > 0 para todo x ∈ σ(A1 ) e y ≤ 0 para todo y ∈ σ(B1 ) , (2.70)
por (2.69) se obtiene que
para todo z ∈ σ −B + α2 A ,

Re(z) > 0 (2.71)
es decir, imponiendo condición (2.70), tenemos garantizada una de las condiciones suficien-
tes para la convergencia de la matriz (2.60). Bajo la condición (2.70), se obtiene
2 mı́n Re(z) : z ∈ σ −B + α2 A

2 Re(z)
mı́n = . (2.72)
z∈σ(−B+α2 A) |z|2 máx {|z|2 : z ∈ σ (−B + α2 A)}
Por (2.71) y (2.72) la condición (2.64) se satisface si
2 mı́n Re(z) : z ∈ σ −B + α2 A

r
< . (2.73)
M2 máx {|z|2 : z ∈ σ (−B + α2 A)}
Puesto que Ak y Bk , k = 1, 2, son matrices hermı́ticas, por el teorema de Gershgorin, véase
el teorema 1.5, si denotamos σ(Ak ) = {λj (Ak ); 1 ≤ j ≤ s}, k = 1, 2, y
Gj (k) = z ∈ C : z − α2 λj (Ak ) ≤ ρ(Bk ) = D α2 λj (Ak ), ρ(Bk ) ,

se sigue que
s
[
2

σ −Bk + α Ak ⊂ Gj (k) , 1 ≤ k ≤ 2. (2.74)
j=1
Por tanto
λmáx −Bk + α2 Ak ≤ α2 λmáx (Ak ) + ρ(Bk ) ,

1 ≤ k ≤ 2, (2.75)
y por (2.67), (2.75) obtenemos
Re(z) ≤ α2 λmáx (A1 ) + ρ(B1 ) , z ∈ σ −B + α2 A ;

(2.76)
Im(z) ≤ α2 λmáx (A2 ) + ρ(B2 ) , z ∈ σ −B + α2 A .

(2.77)
Utilizando en (2.76) y (2.77) que

θ`
0 < α = 2M sen < 2M,
2
se obtiene
2 2
|z|2 ≤ 4M 2 λmáx (A1 ) + ρ(B1 ) + 4M 2 λmáx (A2 ) + ρ(B2 )
. (2.78)
z ∈ σ −B + α2 A
Por (2.64), (2.69), (2.73) y (2.78) se concluye que la matriz dada en (2.60) es convergente
si h i
M 2 4M 2 sen2 θ2` λmı́n (A1 ) − λmáx (B1 )
r< . (2.79)
[4M 2 λmáx (A1 ) + ρ(B1 )]2 + [4M 2 λmáx (A2 ) + ρ(B2 )]2
Resumiendo hemos demostrado el siguiente resultado.
B+B H
TEOREMA 2.2 Sea r > 0 y A, B matrices en Cs×s tales que si B1 = 2 , A1 =
A+AH H A−AH
2 , B2 = B−B
2i , A2 = 2i se cumple
x > 0 , ∀x ∈ σ(A1 ) e y ≤ 0 , ∀y ∈ σ(B1 ) . (2.80)
Entonces para valores suficientemente

pequeños de
r que cumplan la condición (2.79), se
2 θ` 1
obtiene que la matriz I − r 4 sen 2 A − M 2 B es convergente.
Ahora estudiaremos cuándo la expresión matricial dada por (2.59) es convergente, es decir
qué condiciones se han de verificar para que
r
ρ I + 2B < 1 . (2.81)
M
Para que la condición (2.81) se verifique será suficiente que

r
1 + 2 z < 1 , ∀z ∈ σ(B) , (2.82)

M
puesto que por el teorema de la aplicación espectral, véase el teorema 1.6, se sabe que
r n r o
σ I + 2 B = 1 + 2 z / z ∈ σ(B) .
M M
En efecto,
r 2 r r 2
1 + 2 z = 1 + 2 Re(z) + i 2 Im(z)

M M M
r 2 r 2
= 1 + 2 Re(z) + Im(z)
M M2
r 2 2r
= 1 + 4 [Re(z)]2 + [Im(z)]2 + 2 Re(z)
M M
r2 2r
= 1 + 4 |z|2 + 2 Re(z)
M M
r2 2r
≤ 1 + 4 (ρ(B))2 + 2 Re(z) , ∀z ∈ σ(B) .
M M
por tanto, para que se cumpla (2.82) será suficiente imponer las condiciones
2M 2 Re(z)
r<− y Re(z) < 0 , ∀z ∈ σ(B) . (2.83)
(ρ(B))2
Resumiendo, hemos obtenido que las soluciones proporcionadas por el teorema 2.1 son
estables, es decir, permanecen acotadas cuando n → ∞, si los vectores {d` } están acotados,
si se verifican las hipótesis del teorema 2.2 y además se satisfacen las condiciones dadas en
(2.83).
Observación 2.1 El análisis de la estabilidad de las soluciones construidas en el teorema

2.1 se puede hacer, como se verá en el capı́tulo 3 para los sistemas hiperbólicos, utilizando
el concepto de estabilidad en el sentido de estación fija (véase la definición 1.15). Con
este tipo de estabilidad se puede demostrar que las soluciones dadas en el teorema 2.1 son
estables sin imponer restricciones sobre el espectro de las matrices A y B (involucradas en
la solución) ni sobre el parámetro r.
m
En efecto, dado un punto del mallado (X, T ) donde X = M = mh0 , h0 fijo y T = Jk finito,
se debe verificar que las siguientes expresiones estén acotadas:
n n #
I − r 4 sen2 θ2` A − B
M2
, I+ r
M2
B , cuando k −→ 0,
. (2.84)
i.e., n crece, dependiendo del valor de k, con 0 ≤ n ≤ J
Puesto que h0 está fijo y r = hk2 se tiene:

I − r 4 sen2 θ` A − B ≤ 1 + O(k) r

y I + 2 B ≤ 1 + O(k), (2.85)

2 M 2 M
o equivalentemente

I − r 4 sen2 θ` A − B ≤ 1 + S1 k r

y I + 2 B ≤ 1 + S2 k,

2 M2 M
para algunas constantes positivas S1 , S2 .
De las condiciones dadas en (2.85) se obtiene que:
n n #

I − r 4 sen2 θ` A − B
y I + r
B acotadas cuando k −→ 0,
2 M2 M2 ,(2.86)
i.e., n crece, dependiendo del valor de k, con 0 ≤ n ≤ J
ya que
n n
I − r 4 sen2 θ` A − B ≤ I − r 4 sen2 θ` A − B ≤ (1 + O(k))n

2 M2 2 M2
≤ (1 + O(k))J ≤ eJ O(k) ≤ eJck = eT c = O(1),
siendo c una constante positiva; y además
r n r n n
I + 2 B ≤ I + 2 B ≤ (1 + O(k))

M M
≤ (1 + O(k))J ≤ eJ O(k) ≤ eJsk = eT s = O(1),
siendo s una constante positiva.
Hemos visto que las condiciones que aparecen en (2.85) son suficientes para garantizar que
las potencias dadas en (2.84) estén acotadas cuando n crece, a medida que decrece k, hasta
tomar el valor J donde se alcanza el tiempo T , (T = Jk). Por tanto las soluciones dadas
en el teorema 2.1 son estables en el sentido de estación fija.
2.1.5. El problema mixto discreto

En este apartado construiremos una solución exacta del problema (2.7)–(2.10). Observemos
que al separar variables en el problema de contorno (2.7)–(2.9), se llega al problema de
contorno vectorial:
ρ 
−H(m + 1) + 2H(m) − H(m − 1) = − H(m) 
r

H(0) = 0 . (2.87)

CH(M ) + M D [H(M ) − H(M − 1)] = 0

Bajo las hipótesis (2.27) y la existencia de µ ∈ R, µ ∈ σ(−D−1 C), por el teorema de

la aplicación espectral el problema (2.87) es equivalente al siguiente problema de Sturm-
Liouville discreto, véase el apartado 1.2:
ρ

−H(m + 1) + 2H(m) − H(m − 1) = − H(m) , 

r (2.88)
M
H(0) = 0, H(M ) = H(M − 1) , 1 ≤ m ≤ M − 1
M −µ

donde −4r < ρ ≤ 0, entonces por el apartado 2.1.3, las soluciones del problema (2.88) son
H` (m) = sen (mθ` ) d` , d` ∈ Nuc G

e (µ, θ` )
i h
(`−1)π `π
donde θ` ∈ M ,M verifican la ecuación
µ
cos θ` − 1 − M
cot (M θ` ) = . (2.89)
sen θ`
Como los coeficientes del problema (2.88) son escalares, consideramos el problema de Sturm-
Liouville escalar discreto
ρ

−h(m + 1) + 2h(m) − h(m − 1) = − h(m) , 

r (2.90)
M
h(0) = 0 ; h(M ) = h(M − 1) , 1 ≤ m ≤ M − 1 
M −µ

Por el apartado 1.2, el problema (2.90) puede escribirse como un problema de valores
propios matricial
2 −1 0 0 ··· 0
    
h(1) h(1)

 −1 2 −1 0 ··· 0 
 h(2) 


 h(2) 

 0 −1 2 −1 ··· 0  h(3)  ρ h(3) 
.. .. . . .. = −  .
    
.. .. .. ..
. . . r
 

 . . . 
 . 
  . 

 0 0 · · · −1 2 −1  h(M −2)   h(M −2) 
M −2µ
0 0 ··· 0 −1 M −µ h(M −1) h(M −1)
` M −1
El problema (2.90) tiene exactamente M − 1 valores propios reales dados por −ρ
r `=1
,

2 θ` −ρ`
donde ρ` = −4r sen 2 y θ` verifica (2.89). Para cada valor propio r existe una función
propia {h` (m)} = {sen(mθ` )}. Estas funciones propias son ortogonales en
N (1, M − 1) = {k ∈ N : 1 ≤ k ≤ M − 1} ,
con respecto a la función peso 1, véase el apartado 1.2.
Se ha considerado el problema de Sturm-Liouville discreto escalar (2.90), para poder iden-

tificar los vectores d` y d1 que aparecen en las soluciones (2.40) y (2.51) respectivamente.
CASO 1. µ=0 (C singular)

Supongamos que rango G e (0, θ` ) < s , θ` = 2`−1 π , 1 ≤ ` ≤ M − 1. Por superposición
2M −1
de las soluciones del problema de contorno (2.7)–(2.9) obtenidas en el teorema 2.1-(ii), se
sigue que la sucesión vectorial
M −1 n 
θ` B
4 sen2
P
U (m, n) = I −r 2 A− M2
sen (mθ` ) d` , 
`=1
(2.91)



d` e (0, θ` )† G
= I −G e (0, θ` ) S` , S` ∈ Cs ∼ {0}
es una solución del problema (2.7)–(2.9). Imponiendo la condición inicial (2.10) a la solución
(2.91) obtenemos que los vectores d` deben verificar
M
X −1
f (m) = sen (mθ` ) d` , 1 ≤ m ≤ M − 1. (2.92)
`=1
Denotemos por fq (m) la q-ésima componente del vector f (m) ∈ Cs y d`,q la q-ésima com-
ponente del vector d` . De (2.92) se sigue que
M
X −1
fq (m) = sen (mθ` ) d`,q , 1 ≤ m ≤ M − 1. (2.93)
`=1
−1
La expresión (2.93) es una serie de Fourier discreta, y {sen (mθ` )}M `=1 son las funciones
propias del problema de Sturm-Liouville discreto (P.S.L.D.) escalar (2.90). Para determinar
las constantes d`,q , i.e., los coeficientes de Fourier de (2.93), multiplicaremos ambos lados
de (2.93) por sen (mθγ ) para 1 ≤ γ ≤ M − 1, sumaremos el resultado desde m = 1 hasta
m = M − 1 y utilizaremos la ortogonalidad de las funciones propias del P.S.L.D., véase el
apartado 1.2. Es decir:
M −1 M −1 M −1
!
X X X
sen (mθγ ) fq (m) = sen (mθγ ) sen (mθ` ) d`,q
m=1 m=1 `=1
M
X −1

= sen mθγ sen (mθ1 ) d1,q + · · ·
m=1
· · · + sen (mθγ) sen (mθM −1 ) dM −1,q )
M
X −1
= sen2 (mθγ ) dγ,q
m=1
M
X −1
=dγ,q sen2 (mθγ ) , 1 ≤ γ ≤ M − 1,
m=1
y entonces
M −1
P 2`−1
sen m 2M −1 π fq (m)
m=1
d`,q = M −1
2`−1
sen2 m 2M
P
−1 π
m=1
M −1
4 X 2` − 1
= sen m π fq (m) , 1 ≤ q ≤ s ,
2M − 1 2M − 1
m=1
o escrito en forma vectorial

M −1
4 X 2` − 1
d` = sen ν π f (ν) , 1 ≤ ` ≤ M − 1. (2.94)
2M − 1 2M − 1
ν=1
Ya que d` ∈ Nuc G e (0, θ` ) para 1 ≤ ` ≤ M − 1 , en vista de (2.94) se debe verificar la

siguiente condición

2` − 1
f (m) ∈ Nuc G (0, θ` ) ,
e θ` = π, 1 ≤ ` ≤ M − 1. (2.95)
2M − 1
Por tanto bajo la hipótesis (2.95), la expresión
M −1 n 
I − r 4 sen2 θ2` A − B
P
U (m, n) = M2
sen (mθ` ) d` , 
`=1 
, (2.96)
4 M −1  
P 2`−1
d` = sen (νθ` ) f (ν) , θ` = 2M −1 π
2M − 1 ν=1
1 ≤ m ≤ M − 1, n > 0 1 ≤ ` ≤ M − 1
es una solución del problema (2.7)–(2.10). Por la expresión (2.35) la condición (2.95) se
satisface si
f (m) ∈ Nuc G(0) , 1 ≤ m ≤ M − 1, (2.97)

y
Nuc G(0) es subespacio invariante por A + w` B , 1 ≤ ` ≤ M − 1, (2.98)
−1
véase la observación 1.2, donde {w` }M
`=1 se han definido en (2.35). La condición (2.98)
puede escribirse de la forma:

G(0) (A + w` B) I − G(0)† G(0) = 0 , 1 ≤ ` ≤ M − 1, (2.99)
véase la observación 1.3.

CASO 2. µ<0

Sea µ ∈ σ(−D−1 C) y supongamos que rango G(µ, e θ` ) < s donde {θ` }M −1 vienen dadas
`=1
por el teorema 2.1-(iii). Bajo las hipótesis
f (m) ∈ Nuc G(µ) , 1 ≤ m ≤ M − 1, (2.100)
y
Nuc G(µ) es subespacio invariante por A + w` B , 1 ≤ ` ≤ M − 1. (2.101)
−1
donde {w` }M
`=1 se han definido en (2.35), o equivalentemente

G(µ) (A + w` B) I − G(µ)† G(µ) = 0 , 1 ≤ ` ≤ M − 1, (2.102)
entonces una solución del problema (2.7)-(2.10) viene definida por
−1
M n 
θ` B
4 sen2
P
U (m, n) = I −r 2 A− M2
sen (mθ` ) d` , 
`=1 
MP−1 
sen (νθ` ) f (ν) .

(2.103)
d` = ν=1 ,


M −1 
sen2 (νθ` )
P
ν=1
1 ≤ m ≤ M − 1, n > 0, 1 ≤ ` ≤ M − 1
CASO 3. µ=1

Supongamos que rango G e (1, θ` ) < s, 2 ≤ ` ≤ M − 1 y rango R(1)
e < s, entonces por
el teorema 2.1-(iv) la sucesión vectorial
M −1

θ`
n n
A −MB2 r

I − r 4sen2
P
U (m, n)= 2 sen(mθ` )d` + m I + M2 B d1 
`=2 , (2.104)
d` ∈ Nuc G(1,
e θ` ) , d1 ∈ Nuc R(1)
e
1 ≤ m ≤ M − 1, n > 0, 2 ≤ ` ≤ M − 1
es una solución del problema de contorno (2.7)–(2.9). Imponiendo la condición inicial (2.10)
−1
a la solución (2.104) obtenemos que los vectores {d` }M
`=1 deben verificar
M
X −1
f (m) = sen (mθ` ) d` + md1 , 1 ≤ m ≤ M − 1. (2.105)
`=2
Denotemos por fq (m) la q-ésima componente del vector f (m) ∈ Cs y d`,q la q-ésima com-
−1
ponente del vector {d` }M
`=1 . De (2.105) se sigue que
M
X −1
fq (m) = sen (mθ` ) d`,q + md1,q , 1 ≤ m ≤ M − 1. (2.106)
`=2
−1
Para poder identificar los coeficientes {d` }M
`=1 necesitamos tener la ortogonalidad de las
funciones propias
−1
h` (m) = {sen(mθ` )}M
`=2 , (2.107)
y
h1 (m) = m , (2.108)
las cuales son las soluciones de los problemas de Sturm-Liouville (2.90) para los valores
n ρ oM −1 n oM −1
`
propios − = 4 sen2 θ2` y 0 respectivamente. Observemos que los M − 1
r `=2 `=2
valores propios son distintos entre sı́.
Por la teorı́a de Sturm-Liouville, véase el apartado 1.2, sabemos que las funciones propias
dadas en (2.107) son ortogonales en N (1, M − 1) respecto a la función peso 1, por tanto
para garantizar que el conjunto de las funciones propias:

m, para ` = 1
h` (m) =
sen(mθ` ) , para 2 ≤ ` ≤ M − 1
es ortogonal en N (1, M − 1) respecto a la función peso 1, sólo tendremos que estudiar

si una función cualquiera de (2.107) y la función definida en (2.108) son ortogonales en
N (1, M − 1) a la función peso 1.
En efecto, puesto que una función h` (m) de (2.107) y h1 (m) son funciones propias del
problema (2.90) verifican su ecuación
ρ
− sen ((m + 1)θ` ) + 2 sen (mθ` ) − sen ((m − 1)θ` ) + sen (mθ` ) = 0 , (2.109)
r
−(m + 1) + 2m − (m − 1) = 0 , (2.110)
Multiplicando la ecuación (2.109) por h1 (m) y (2.110) por h` (m), y restando las ecuaciones
resultantes obtenemos:
ρ
m sen (mθ` ) = −m [sen ((m + 1)θ` ) + sen ((m − 1)θ` ) − 2 sen (mθ` )] .
r
1≤m≤M −1
Si sumamos m, desde 1 hasta M − 1, para cada miembro de la igualdad anterior obtenemos

M −1 M −1
ρ X X
m sen (mθ` ) = sen (mθ` ) (−2m + 2m) = 0
r
m=1 m=1
ρ
y puesto que r 6= 0 se obtiene
M
X −1
m sen (mθ` ) = 0
m=1
es decir, son ortogonales en N (1, M − 1) respecto a la función peso 1. Ahora ya podemos

−1
pasar a determinar los vectores {d` }M
`=1 . Para ello procederemos de la misma forma que
−1
se vio en el Caso 1. Por un lado obtenemos la expresión de los vectores {d` }M
`=2
M −1 M −1 M −1 M −1
!
X X X X
sen (mθγ ) fq (m) = sen (mθγ ) sen (mθ` ) d`,q +d1,q m sen (mθγ )
m=1 m=1 `=2 m=1
M
X −1
= dγ,q sen2 (mθγ ) , 1 ≤ γ ≤ M − 1,
m=1
es decir
MP−1
sen (mθ` ) fq (m)
m=1
d`,q = M −1
, 1 ≤ q ≤ s,
sen2 (mθ` )
P
m=1

MP−1
sen (νθ` ) f (ν)
ν=1
d` = M −1
, 2 ≤ ` ≤ M − 1. (2.111)
sen2 (νθ
P
`)
ν=1
Por otro lado, para obtener el vector d1 procedemos de la siguiente forma

−1 −1 M −1
M M
! M −1
X X X X
m fq (m) = m sen (mθ` ) d`,q + m2 d1,q
m=1 m=1 `=2 m=1
M −1
M3 M2

X
2 M
= d1,q m = − + d1,q .
3 2 6
m=1
Por tanto
M −1
1 X
d1,q = m fq (m) 1 ≤ q ≤ s,
M3 M2 M
3 − 2 + 6 m=1

M −1
1 X
d1 = ν f (ν) . (2.112)
M3 M2 M
3 − 2 + 6 ν=1
Además como se tiene que verificar que d` ∈ Nuc G e (1, θ` ), 2 ≤ ` ≤ M − 1 y que d1 ∈

Nuc R(1),
e en vista de (2.111) y (2.112) será suficiente exigir las siguientes condiciones
f (m) ∈ Nuc G(1) , (2.113)

Nuc G(1) subespacio invariante por las matrices A + w` B y B , (2.114)

M −1 M −1
donde {w h definidos en (2.35) y {θ` }`=2 son las soluciones de la ecuación (2.89)
i ` }`=2 están
en J` = (`−1) `
M π, M π para 2 ≤ ` ≤ M − 1.
Las condiciones dadas en (2.114) son equivalentes a

G(1) (A + w` B) I − G(1)† G(1) = 0 , 2 ≤ ` ≤ M − 1, (2.115)
y
G(1)B I − G(1)† G(1) = 0 , (2.116)
respectivamente.
CASO 4. µ ∈ ]0, +∞[, µ 6= 1
En este caso habı́amos probado analı́ticamente la existencia de M −2 soluciones, θ` ∈ J` ,

2 ≤ ` ≤ M − 1, para la ecuación (2.89) y por el teorema 2.1-(v) obtuvimos que bajo la
e (µ, θ` ) < s, la sucesión {U` (m, n)}M −1 dada por (2.40) definı́a
hipótesis de que rango G `=2
M −2 soluciones no triviales del problema (2.7)–(2.9). Ahora bien, notemos que la teorı́a de
Sturm-Liouville para problemas de valores propios
matriciales nos garantiza la existencia
ρ`
de M − 1 valores propios, − r = −4 sen θ2` en ]0, π[. Y puesto que las ecuaciones
2
(2.89) y h(M ) = MM−µ h(M −1) siendo h(m) = sen(mθ` ) son la misma, se garantiza que en
π
el primer subintervalo, ]0, M [, también hay solución para la ecuación (2.89) y por tanto
existen M − 1 soluciones para el problema de contorno (2.7)–(2.9). Por superposición de
las M −1 soluciones y bajo las hipótesis

rango G(µ,
e θ` ) < s , 




f (m) ∈ Nuc G(µ) ,  1 ≤ ` ≤ M − 1 , 1 ≤ m ≤ M − 1,




†

G(µ)(A + w` B) I − G(µ) G(µ) = 0

obtenemos la siguiente solución para el problema mixto (2.7)–(2.10)
−1
M n 
θ` B
4 sen2
P
U (m, n) = I −r 2 A− M2
sen (mθ` ) d` , 
`=1 
MP−1 
sen (νθ` ) f (ν) .

d` = ν=1 ,


M −1 
sen2 (νθ` )
P
ν=1
1 ≤ m ≤ M − 1, n > 0, 1 ≤ ` ≤ M − 1
Notemos que las condiciones que se han exigido para la existencia de solución del problema
mixto (2.7)–(2.10) son las mismas que para el caso µ < 0, al igual que sucede con la solución
construida.
Observación 2.2 A la vista de los resultados obtenidos para cualquier valor propio real,
µ, de la matriz −D−1 C, podemos establecer una cota más concreta del parámetro r que la
dada en (2.79). Utilizando que sen θ es creciente en 0, π2 concluimos que

2 2 θ` 2 2 θ1
4M sen λmı́n (A1 ) − λmáx (B1 ) ≥ 4M sen λmı́n (A1 ) − λmáx (B1 )
2 2
(si µ 6= 1)

2 2 θ` 2 2 θ2
4M sen λmı́n (A1 ) − λmáx (B1 ) ≥ 4M sen λmı́n (A1 ) − λmáx (B1 )
2 2
(si µ = 1)
Por tanto obtenemos las siguientes cotas de r:
• Si µ 6= 1 h i
M 2 4M 2 sen2 θ21 λmı́n (A1 ) − λmáx (B1 )
r< . (2.117)
• Si µ = 1 h i
M 2 4M 2 sen2 θ22 λmı́n (A1 ) − λmáx (B1 )
r< . (2.118)

TEOREMA 2.3 Con la notación y las hipótesis de los teoremas 2.1 y 2.2, asumamos las
condiciones dadas en (2.80) y que {f (m)} está acotado.

−1
(i) Si C es singular, µ = 0, rango G e 0, 2`−1 π
2M −1 < s, 1 ≤ ` ≤ M −1, y {f (m)}Mm=1
y las matrices A + w` B verifican las condiciones (2.97) y (2.99) respectivamente,
entonces si r verifica (2.117), {U (m, n)} definida por (2.96) es una solución exacta
estable del problema mixto (2.7)–(2.10).

(ii) Si µ = 1, rango G e (1, θ` ) < s, donde {θ` }M −1 son las soluciones de (2.89) en
`=2
i h
(`−1)π ` −1
J` = M , M π , 2 ≤ ` ≤ M − 1, rango R(1)
e < s, {f (m)}M
m=1 verifica (2.113), y
−1
las matrices {A + w` B}M`=2 y B verifican las condiciones (2.115) y (2.116) respec-
tivamente, entonces si además se verifica la condición (2.83) y r satisface (2.118),
obtenemos que {U (m, n)} definida por (2.104) es una solución exacta estable del pro-
−1
blema mixto (2.7)–(2.10), siendo {d` }M `=2 y d1 los vectores definidos en (2.111) y
(2.112) respectivamente.

(iii) Si µ ∈ R ∼ {0, 1}, rango G e (µ, θ` ) < s, donde {θ` }M −1 son las soluciones de (2.89)
`=1
i h
(`−1)π ` M −1 M −1
en J` = M , M π , 1 ≤ ` ≤ M − 1, y {f (m)}m=1 y las matrices {A + w` B}`=1
verifican las condiciones (2.100) y (2.102) respectivamente, entonces si r verifica
(2.118), {U (m, n)} definida por (2.103) es una solución exacta estable del problema
mixto (2.7)–(2.10).
2.1.6. El método de las proyecciones

En este apartado abordaremos la construcción de soluciones del problema (2.7)–(2.10) para
una función f (m) satisfaciendo condiciones más generales que las exigidas en el apartado
2.1.5. Supongamos que
{µ1 , . . . , µq } ⊂ σ −D−1 C ∩ R ,

(2.119)
donde µj 6= µk para 1 ≤ j, k ≤ q, j 6= k. Introducimos la matriz G(µj ) en Cs×s
G(µj ) = D−1 C + µj I , 1 ≤ j ≤ q. (2.120)
Como los polinomios x − µj son mútuamente coprimos (primos 2 a 2), por el teorema de
descomposición, véase el teorema 1.9, si L(x) es el polinomio
L(x) = (x − µj ) (x − µ1 ) (x − µ2 ) · · · (x − µj−1 ) (x − µj+1 ) · · · (x − µq ) , (2.121)
entonces
S = Nuc L −D−1 C = Nuc G (µ1 ) ⊕ Nuc G (µ2 ) ⊕ · · · ⊕ Nuc G (µq ) .

(2.122)
Consideremos la sucesión de polinomios coprimos de grado q − 1 definidos por

q
Y
Qj (x) = (x − µk ) , 1≤j≤q. (2.123)
k=1,k6=j
Entonces por el teorema de Bezout, véase el teorema 1.8, tomando los escalares αj de la
forma
 −1
q
 Y 
αj = (µj − µk ) , 1≤j≤q,
 
k=1,k6=j
obtenemos
q
X
1 = Q(x) = αj Qj (x) . (2.124)
j=1
Notemos que Q(x) es el polinomio de interpolación de Lagrange satisfaciendo Q(µj ) = 1

para 1 ≤ j ≤ q.
Considerando el polinomio
T (µj ) = G (µ1 ) · · · G (µj−1 ) G (µj+1 ) · · · G (µq ) , (2.125)
y aplicando el cálculo funcional matricial sobre la matriz −D−1 C, de (2.120), (2.121),

(2.124) y (2.125) se sigue que
I = Q −D−1 C

Xq
αj Qj −D−1 C

=
j=1
q
X
q−1
= (−1) αj G(µ1 )G(µ2 ) · · · G(µj−1 )G(µj+1 ) · · · G(µq )
j=1
Xq
= (−1)q−1 αj T (µj ) . (2.126)
j=1
Dado el problema (2.7)–(2.10) consideramos la proyección de la sucesión {f (m)}M

m=0 sobre
el subespacio Nuc G(µj ) definida por
fbj (m) = (−1)q−1 αj T (µj )f (m) , 1 ≤ j ≤ q, 0 ≤ m ≤ M , (2.127)
y observar que por (2.126) y (2.127) se sigue que
q
X q
X
fbj (m) = (−1)q−1 αj T (µj )f (m)
j=1 j=1
 
q
X
= (−1)q−1 αj T (µj ) f (m) (2.128)
j=1
= If (m) = f (m) . (2.129)
Supongamos que f (m) ∈ S, 0 ≤ m ≤ M , i.e.,
L(−D−1 C)f (m) = 0 , 0≤m≤M, (2.130)
entonces por (2.121), (2.127) y (2.130) obtenemos
G(µj )fbj (m) = (−1)q−1 αj G(µj )T (µj )f (m)

= (−1)q−1 (−1)q αj L(−D−1 C)f (m) = 0
1 ≤ j ≤ q, 0 ≤ m ≤ M ,
o equivalentemente
fbj (m) ∈ Nuc G(µj ) , 1 ≤ j ≤ q, 0 ≤ m ≤ M . (2.131)
Por tanto se consideran q problemas diferentes obtenidos al reemplazar la condición inicial

(2.10) por
U (m, 0) = fbj (m) , 0 ≤ m ≤ M , 1 ≤ j ≤ q. (2.132)
Sea (Pj ), 1 ≤ j ≤ q, el problema definido por (2.7)–(2.9) y (2.132). Observar que por resul-
tados previos, el cambio de la condición inicial en cada problema (Pj ) sólo modifica a los
vectores {d` } que aparecen en (2.94), (2.103), (2.111) ón (2.112)

o porque ahora los coeficien-
tes de Fourier discretos involucran a las proyecciones fbj (m) en lugar de al vector de la
condición inicial {f (m)}. Notemos que para los problemas (Pj ), los coeficientes de Fourier
toman la forma
si µj ∈ R ∼ {1}:
M −1
P (j)
sen νθ` f (ν)
(j) ν=1
d` = (−1)q−1 αj T (µj ) M −1 , 1 ≤ ` ≤ M − 1, (2.133)
(j)
sen2 νθ`
P
ν=1
o si µj = 1
M −1 
P (j)
sen νθ` f (ν)



(j)

ν=1
(−1)q−1 αj ≤ ≤ −

d` = T (µj ) M −1 , 2 ` M 1


(j)

sen2 νθ`
P 
ν=1
(2.134)




M −1


d1 = (−1)q−1 αj T (1) M 3
(j) 1 P 
νf (ν)


2
− M2 + M 
3 6 ν=1
donde el superı́ndice (j) denota que se está trabajando con el valor propio µj . Si denotamos
M −1
P (j)
sen νθ` f (ν)
(j) ν=1
d` (f (ν)) = M −1 ,
(j)
sen2
P
νθ`
ν=1
M −1
(j) 1 X
d1 (f (ν)) = νf (ν),
M3 M2 M
3 − 2 + 6 ν=1
n o
entonces los coeficientes de Fourier que involucran a fbj (m) dados por (2.133) y (2.134)
pueden escribirse en función de {f (m)} por medio de la expresión
si µj ∈ R ∼ {1}:
(j) (j)
d` = (−1)q−1 αj T (µj ) d` (f (ν)) , 1 ≤ ` ≤ M − 1, (2.135)
o si µj = 1:

(j) (j)
d` = (−1)q−1 αj T (1)d` (f (ν)) , 2 ≤ ` ≤ M − 1 

. (2.136)
(j) (j)
(−1)q−1 αj T (1)d1 (f (ν))

d1 = 

TEOREMA 2.4 Con la notación del teorema 2.3, sean µ1 , . . . , µq valores propios reales
y distintos de la matriz −D−1 C. Sean G(µj ), L(x) y T (µj )matrices
en C
s×s definidas por
(2.120), (2.121) y (2.125) respectivamente. Sean R(1)

e yG e µj , θ(j) las matrices definidas
`
por
 

G(1)
 G(µj )
(j)
G(µj ) A + w` B
 
G(1)B
   
e µj , θ(j)
 
R(1)
e =

.. ,

G ` = .. ,
 .  
 . 

p(j) −1
G(1)B t−1
 
(j) `
G(µj ) A + w` B
(j) (j)
donde t y p` es el grado del polinomio minimal de las matrices B y An + w` oB, respec-
tivamente. Supongamos que {f (m)} es acotado y verifica (2.130) (i.e. fbj (m) verifican
(2.131)), y que las condiciones dadas en (2.80) se satisfacen.
(j) 1
(i) Si µj ∈ R ∼ {1}, definimos w` = − (j) , 1 ≤ ` ≤ M − 1, donde
θ
4M 2 sen2 `2

(j) 2` − 1 (j)
i h
θ` = π si µj = 0, y θ` ∈ (`−1)π
M , `π
M , 1 ≤ ` ≤ M −1, son las soluciones
2M − 1
de la ecuación

cos θ`(j) − 1 − µMj

(j)
cot M θ` = ,1 ≤ ` ≤ M − 1.
(j)
sen θ`

En este caso, supongamos que rango G e µj , θ(j) < s, 1 ≤ ` ≤ M−1, que las matrices
`
(j)
A + w` B verifican

(j)
G(µj ) A + w` I − G (µj )† G (µj ) = 0 , 1 ≤ ` ≤ M − 1.
y que el parámetro r satisface (2.117).

(j) 1
(ii) Si µj = 1, definimos w` = − (j)
, 2 ≤ ` ≤ M − 1, donde
θ`
4M 2 sen2 2
i h
(j) (`−1)π `π
θ` ∈ M ,M , 2 ≤ ` ≤ M − 1, son las soluciones de la ecuación

cos θ`(j) − 1 − µMj

(j)
cot M θ` = , 2 ≤ ` ≤ M − 1.
(j)
sen θ`

En este caso, supongamos que rango G e 1, θ(j) < s, para 2 ≤ ` ≤ M − 1 y que
`
n oM −1
(j)
rango R(1)
e < s; que se verifica la condición (2.83); las matrices A + w` B
`=2
y B verifican las condiciones

(j)
G(1) A + w` B I − G(1)† G(1) = 0 ,

2 ≤ ` ≤ M − 1,
G(1)B I + G(1)† G(1) = 0 ;

y el parámetro r satisface (2.118).
Entonces bajo las hipótesis consideradas anteriormente en (i)-(ii) se obtiene que

q
X
U (m, n) = Uj (m, n) , 1 ≤ m ≤ M − 1, n > 0, (2.137)
j=1
es una solución convergente (por ser consistente y estable) a la solución del problema mixto
continuo (2.1)–(2.4), donde
si µj ∈ R ∼ {1}:
M −1
" (j)
! !#n
X
2 θ` B
(j)

(j)
Uj (m, n) = I − r 4 sen A− 2 sen m θ` d` ,
2 M
`=1
(j) M −1
n o
siendo d` los vectores definidos en (2.135);
`=1
o si µj = 1:
M −1 (j)
! !!n
θ` B r n
(j) (j) (j)
X
Uj (m, n)= I − r 4sen2 A− 2 sen mθ` d` + m I 2 B d1
2 M M
`=2
(j) M −1
n o
(j)
siendo d` y d1 los vectores definidos en (2.136).
`=2
2.1.7. Ejemplos
A continuación presentamos dos ejemplos para ilustrar el funcionamiento del método y un
tercero, a modo de problema test, para comparar la solución teórica exacta con la solución
numérica construida.
EJEMPLO 2.1 Consideremos el problema (2.1)-(2.4) con las matrices:

 11
− 25 3
 13
− 6 − 61
   
2 2 2 0 −1 0 0
A =  −9 79 −9 , B =  1 −1 0 , C =  4 −1 −2 ,
1 1
−3 25 1 3 3 −2 3 − 12 −1
 
1 0 0
D =  0 2 1  , F (mh) = f (m) = (f1 (m), f2 (m), f3 (m))T ∈ C3
0 1 1
1
h= M , 1 ≤ m ≤ M − 1.
Tenemos  
1 0 0
−D−1 C =  −1 21 1 ,
−2 0 0
por tanto σ(−D−1 C) = 0, 12 , 1 . Denotamos µ1 = 0, µ2 = 12 , µ3 = 1.

Veamos si se cumplen las siguientes condiciones:

(a) rango G e µj , θ(j) < 3, para cada 1 ≤ j ≤ 2 con 1 ≤ ` ≤ M − 1 ;
`
(b)



 rango R(1)
e < 3;

e 1, θ(3)

 rango G < 3, 2 ≤ ` ≤ M − 1 .

`
En nuestro caso tenemos
 (j) 
11 13 (j) w` 25 3
2 − 6 w` − 6 − 2 2
 
 
(j)  (j) (j)

,1 ≤ j ≤ 3 .
A + w` B = 
 w` − 9 79 − w` −9 
 
 
(j) (j)
w` w` (j)
3 −3 3 + 25 1 − 2w`
(j) (j)
El polinomio minimal de las matrices A+w` B , ∀w` , 1 ≤ j ≤ 3, coincide con su polinomio
(j)
caracterı́stico, por tanto, p` = 3, ∀j tal que 1 ≤ j ≤ 3. Por otro lado tenemos que el
polinomio minimal de B es de grado 2, es decir, t = 2. Por tanto tendremos que estudiar
el rango de las siguientes matrices:
 
G(µj )
 
 
 G(µ ) A + w(j) B
e µj , θ(j)

G =
 j `  , 1 ≤ j ≤ 3,
` 
 
 2 
(j)
G(µj ) A + w` B
y
 
G(1)
R(1)
e = .
G(1)B
1
Se comprueba que los valores propios µ1 = 0 y µ = 2 no verifican la condición (a).
Para µ3 = 1 tenemos:  
0 0 0
G(1) = D−1 C + I =  1 1
2 −1
,
2 0 1
 
0 0 0
 
 
1

 1 2 −1 

 
 

 2 0 1 

 
 

 0 0 0 

 
 
G(1,
e θ` ) =  2 (2 − w` ) 2 − w` 2w` − 4 ,


 
 
8 − 4w` 0 4 − 2w`
 
 
 
 
 

 0 0 0 

 
 4 w2 − 4w` + 4 2 w2 − 4w` + 4 −4 w2 − 4w` + 4
 

 ` ` ` 
 
8 w`2 − 4w` + 4 4 w`2 − 4w` + 4

0
 
0 0 0
 
 
1

 1 2 −1 

 
 
2 0 1
 
 
R(1) = .
e  
0 0 0
 
 
 
 
−2 −1
 
 2 
 
 
−4 0 −2

e (1, θ` ) = 2 < 3, ∀w` , 2 ≤ ` ≤ M − 1 y rango R(1)
Por tanto, rango G e = 2 < 3. Además
 12 14 
0 41 41
 
† 10 2
 
 0
G(1) =  41 − 41 ,

 
0 − 24
41
13
41
y Nuc G(1) es subespacio invariante por las matrices A + w` B, 2 ≤ ` ≤ M − 1, y B, es decir

G(1) (A + w` B) I − G(1)† G(1) = 0 , 2 ≤ ` ≤ M − 1 ,
y
G(1)B I − G(1)† G(1) = 0 .
Para que el vector {f (m)} ∈ Nuc G(1) debemos exigir
 
f1 (m)
f (m) =  −6f1 (m)  ,
−2f1 (m)
y por el teorema 2.3-(ii) la solución numérica del problema viene dada por
M −1 n
X
2 θ` B r n
U (m, n)= I − r 4sen A− 2 sen(mθ` )d` + m I + 2 B d1
2 M M
`=2
donde
MP−1
sen (νθ` ) f (ν)
ν=1
d` = M −1
, 2 ≤ ` ≤ M − 1.
sen2 (νθ
P
`)
ν=1
i h
(`−1)π `π
con θ` ∈ M ,M , 2 ≤ ` ≤ M − 1 , la solución de
1

cos (θ` ) − 1 − M
cot (M θ` ) = ,
sen (θ` )
y
M −1
1 X
d1 = νf (ν)
M3 M2 M
3 − 2 + 6 ν=1
En este caso tenemos que las matrices

 11
− 43 − 34

2 4
 
 43 
A1 =  − 4
 79 8   ; σ(A1 ) = {4, 81,3417, 0,1583} ,
 
3
−4 8 1
− 13 5 1
 
6 12 6
 
5 1
 
B1 = 
 12 −1 6
 ; σ(B1 ) = {−2, −0,8287, −2,3379} ,

 
1 1
6 6 −2
y además
7
σ(B) = −2, − ,
6
es decir, se verifican las condiciones dadas en (2.80) y (2.83). Entonces tomando {f (m)}
acotado y r satisfaciendo simultáneamente las condiciones (2.83) y (2.118) se garantiza la
estabilidad de la solución construida.
EJEMPLO 2.2 Consideremos el problema (2.1)–(2.4) con los datos
− 21
 
    0 0
2 0 0 −2 0 0 
 1


A =  1 1 0 , B =  −1 −1  −2
0 , C =  0 0 
,
0 0 2 0 0 −5  
1
2 0 3
 
1 0 0
D =  2 −1 0  , F (mh) = f (m) = (f1 (m), f2 (m), f3 (m))T ∈ C3 ,
−1 0 3
1
h= M , 1 ≤ m ≤ M − 1.
En este caso tenemos

 1 
2 0 0
 
−1 1
 
−D C=
 2 0 0 ;
 
0 0 −1

−1 1 1
σ(−D C) = −1, 0, , µ1 = 0, µ2 = , µ3 = −1 ,
2 2
y con la notación utilizada en la teorı́a tomando µ1 y µ2 obtenemos
− 21
   
0 0 0 0 0
   
 1   1 1 
 −2 0 0
G (µ1 ) =   −2 2 0 
; G (µ2 ) = 
 .
   
3
0 0 1 0 0 2
Notemos que para cada µj tenemos
 
G(µj )
 
 
 G(µ ) A + w(j) B
e µj , θ(j)

G =
 j `  ,1 ≤ j ≤ 2 , 1 ≤ ` ≤ M −1 .
` 
 
 2 
(j)
G(µj ) A + w` B
Para µ1 = 0 se tiene
− 12
 
0 0
 
 

 − 12 0 0 

 
 

 0 0 1 

 
 
(1)
w` − 1

0 0

 
 
 
(1)
 
(1) w` − 1 0 0
G
e µ1 , θ =  , 1≤`≤M −1.
 
`  
 
(1)
0 0 2 − 5w`
 
 
 
 
 −2(w(1) − 1)2
 
 ` 0 0 

 
 
(1)
 −2(w` − 1)2 0 0
 

 
 
(1)
0 0 (5w` − 2)2

e µ1 , θ(1) = 2 < 3 para cualquier w(1) con 1 ≤ ` ≤ M −1.
Por tanto, rango G ` `
Además
 
−1 −1 0
G(µ1 )† =  0 0 0  ,
0 0 1
(1)
y Nuc G(µ1 ) es un subespacio invariante por la matriz A + w` B , para 1 ≤ ` ≤ M − 1 ,
porque

(1)
G(µ1 ) A + w` B I − G(µ1 )† G(µ1 ) = 0 .
1
Para µ2 = 2 se tiene
 
0 0 0
 
 

 − 21 1
2 0 

 
 
3

 0 0 2


 
 

 0 0 0 

 
 
 (2) (2)
1−w` 1−w`

 − 2 2 0 
e µ2 , θ(2) = 
 
G  , 1≤`≤M −1.
` 
 (2)


 3 2−5w` 
 0 0 2

 
 
 

 0 0 0 

 2 2
 

(2) (2)
 w` −1 w` −1 
 − 0
 
 2 2 

(2) 2 
 

3 5w` −2
0 0 2
Para este valor propio obtenemos
 
0 −1 0
G(µ2 )† =  0 1 0  ,
0 0 23
(2)
y que Nuc G(µ2 ) es un subespacio invariante por A + w` B , 1 ≤ ` ≤ M − 1 , porque

(2)
G(µ2 ) A + w` B I − G(µ2 )† G(µ2 ) = 0 .
Consideremos el teorema 2.4 con µ1 = 0, µ2 = 12 , q = 2 . Para garantizar que f (m) ∈ S =

Nuc L(−D−1 C) = Nuc G(µ1 ) ⊕ Nuc G(µ2 ), imponemos la condición equivalente de que las
proyecciones fb1 (m) ∈ Nuc G(µ1 ) y fb2 (m) ∈ Nuc G(µ2 ). Teniendo en cuenta las expresiones
de G(µ1 ) y G(µ2 ) es fácil comprobar que
   b 
0 f2,1 (m)
   
   
fb1 (m) =  f1,2 (m) 
b  y  fb2,1 (m)  .
fb2 (m) =  
   
0 0
Por tanto  
fb2,1 (m)
 
 
f (m) = fb1 (m) + fb2 (m) = 
 fb1,2 (m) + fb2,1 (m)
,

 
0
y por el teorema 2.4 la solución numérica del problema viene dada por
2
X
U (m, n) = Uj (m, n) ,
j=1
donde
M −1 h in
2`−1 B 2`−1 (1)
I −r 4A sen2 2(2M
P
U1(m, n) = −1) π − M2
sen m 2M −1 π d`
`=1
1 ≤ m ≤ M −1 , n > 0,
 
0 0 0
  M −1
(1) 8  1
 − 1
X 2` − 1
d` = 2 2 0  sen ν π f (ν) ,
2M − 1 


 ν=1 2M − 1
3
0 0 2
1 ≤ ` ≤ M − 1,
M −1
(2) n
P 2 θ` B (2) (2)
U2 (m, n) = I − r 4 A sen 2 − M2
sen mθ` d` ,
`=1
1 ≤ m ≤ M − 1, n > 0,
i h
(2) (`−1)π `π
con θ` ∈ M ,M , 1 ≤ ` ≤ M − 1 , la solución de

(2) µ2
cos θ` − 1− M
(2)
cot M θ` = ,
(2)
sen θ`
y
 MP−1
(2)
sen νθ` f (ν)

1 0 0
(2) ν=1
d` = 1 0
 0 
M −1 , 1≤`≤M −1.
0 0 −2 (2)
sen2 νθ
P
`
ν=1
En este caso tenemos

 1 
2 2 0
 
 1
 2174 781
A1 =  2 1 0 
; σ(A1 ) = , ,2 ,

  985 985
0 0 2
−2 − 21
 
0
 
 1  2174 781
 −2
B1 =  −1 0 
; σ(B1 ) = −
985
,−
985
, −5 ,
 
0 0 −5
y tomando {f (m)} acotado y r satisfaciendo (2.117) se garantiza la estabilidad de la solu-
ción construida.
Observación 2.3 Es fácil comprobar que para el valor propio µ3 = −1 no se satisfacen

las hipótesis (i) del teorema 2.4.
EJEMPLO 2.3 Consideremos el problema (2.1)-(2.4) con las matrices:

     
2 0 0 −1 0 0 −1 0 0
A =  0 1 1 , B =  0 −1 0 , C =  0 −1 0 ,
0 0 1 0 0 −1 0 0 −1
 
1 0 0
D =  0 1 0 , y el vector F (x) = (x, x, x)T ∈ C3 .
0 0 1
Fácilmente se prueba que
u(x, t) = (xe−t , xe−t , xe−t )T , (2.138)
es la única solución teórica del problema dado.
Ahora vamos a utilizar los resultados teóricos obtenidos para calcular una solución numéri-
ca del problema dado. Posteriormente compararemos ambas soluciones.
Tenemos:
 
1 0 0
−D−1 C =  0 1 0 ,
0 0 1
por tanto σ(−D−1 C) = {1}. Denotamos µ = 1. Estamos en el caso del teorema 2.3-(ii).
Veamos que se cumplen las siguientes condiciones:
(a)



 rango G
e (1, θ` ) < 3, 2 ≤ ` ≤ M − 1;

 rango R(1)
 e < 3;
(b)
 G(1) (A + w` B) I − G(1)† G(1) = 0,

2 ≤ ` ≤ M − 1;
G(1)B I − G(1)† G(1) = 0 .


En nuestro caso tenemos
 
2 − w` 0 0
 
 
A + w` B = 
 0 1 − w` 1 , 2 ≤ ` ≤ M − 1 .

 
0 0 1 − w`
El polinomio minimal de las matrices A + w` B coincide con su polinomio caracterı́stico,

por tanto, p` = 3, para 2 ≤ ` ≤ M − 1. Por otro lado tenemos que el polinomio minimal
de B es de grado 1, i.e., t = 1. Por tanto tendremos que estudiar el rango de las siguientes
matrices:
 
G(1)
 
 
 , 2 ≤ ` ≤ M −1 ,
 G(1) (A + w` B)
G (1, θ` ) = 
e

 
G(1) (A + w` B)2
y
 
0 0 0
R(1)
e = G(1) = D−1 C + I =  0 0 0  .
0 0 0

Obtenemos que rango G e (1, θ` ) = rango R(1)
e = 0 < 3. Además puesto que G(1) es la
matriz nula, se verifican las condiciones de (b). Notemos que el vector de la condición inicial
f (m) = (mh, mh, mh)T está acotado y pertenece al Nuc G(1). Por tanto, por el teorema
2.3-(ii) la solución numérica del problema viene dada por
M −1 n
X
2 θ` B r n
U (m, n)= I − r 4sen A− 2 sen(mθ` )d` + m I + 2 B d1
2 M M
`=2
1 ≤ m ≤ M − 1, n > 0
donde
MP−1
sen (νθ` ) f (ν)
ν=1
d` = M −1
, 2 ≤ ` ≤ M − 1,
sen2 (νθ` )
P
ν=1
i h
(`−1)π `π
siendo θ` ∈ M ,M , 2 ≤ ` ≤ M − 1 , la solución de
1

cos (θ` ) − 1 − M
cot (M θ` ) = , (2.139)
sen (θ` )
y
M −1
1 X
d1 = νf (ν) .
M3 M2 M
3 − 2 + 6 ν=1
Notemos que
 
2 0 0
 
1  ; σ(A1 ) = 2, 1 , 3 ,
 
 0 1
A1 =  2  2 2
 
0 21 1
 
−1 0 0
 
 
B1 = 
 0 −1  ; σ(B1 ) = {−1} ; σ(B) = {−1} .
0 
 
0 0 −1
Para que la solución numérica construida sea estable se han de verificar simultáneamente
las siguientes condiciones sobre el parámetro r > 0:
h i
M 2 −λmáx (B1 ) + 4M 2 sen2 θ22 λmı́n (A1 )
r< , (2.140)
π 2π

donde θ2 es la única solución de la ecuación (2.139) en M, M , que es el segundo subin-
tervalo de ]0, π[, y
2M 2 Re(z)
r<− , ∀z ∈ σ(B) . (2.141)
(ρ(B))2
Para poder comparar numéricamente la solución construida con la solución teórica exacta,
tomaremos distintos valores de M .
Caso 1. M = 10

θ2
π 2π
Para la elección de r en (2.140) necesitamos averiguar el valor de sen2
2 , θ2 ∈ M , M .
M −1
Sabemos que − ρr` `=2 = 4 sen2 θ2` son los valores propios del problema escalar de

Sturm-Liouville (2.90) para µ = 1, i.e., son los valores propios de la matriz tridiagonal
e (M − 1 × M − 1), de la forma:
A,
2 −1 0 0 ··· 0
 

 −1 2 −1 0 ··· 0 

 0 −1 2 −1 ··· 0 
A
e = .. .. . . .. , (9 × 9),

.. ..
. . .


 . . . 

 0 0 · · · −1 2 −1 
M −2µ
0 0 ··· 0 −1 M −µ
que se pueden
calcular de forma eficiente con elpaquete informático
Matlab.
Por otro
θ` θ`
lado, sen 2 es creciente en ]0, π[ =⇒ sen 2 ≥ sen 2 > sen θ21 , ∀` tal que
2 2 2 θ2 2

2 ≤ ` ≤ M − 1 =⇒ 4 sen2 θ22 es el segundo valor propio más pequeño de la matriz A, e es
decir,

2 θ2
4 sen = 0,2199803986793. (2.142)
2
Sustituyendo (2.142) en la expresión (2.140) y teniendo en cuenta (2.141) se obtiene que r
debe satisfacer simultáneamente las condiciones:
r < 0,00176041700848 y r < 200 ,

donde se ha tenido en cuenta que
 
0 0 0
A − AH
 
 1
 1 1
A2 =  0
= 0 − 2i ; σ(A2 ) = − , 0, ,
2i 

 2 2
1
0 2i 0
B − BH
B2 = = 0; σ(B2 ) = {0} .
2i
Tomemos r = 0,0017. Por tanto, el paso temporal k será k = 0,000017. Para calcular
−1
numéricamente la solución U (m, n) necesitamos conocer los valores de {θ` }M `=2 , ya que
−1
en U (m, n) nos aparece el término sen (mθ` ). Por teorı́a sabemos que la sucesión {θ` }M`=2
verifica la ecuación
2r + ρ`
cos (θ` ) = , 2 ≤ ` ≤ M − 1,
2r
de donde obtenemos que
ρ
1 `
θ` = arc cos − −2 + − , 2≤`≤M −1, (2.143)
2 r
M −1
siendo − ρr` `=2 los valores propios de la matriz A.

e Se comprueba fácilmente que los θ`
obtenidos utilizando la ecuación (2.143) son las soluciones de la ecuación (2.139). Notar
−1
que se ha utilizado este procedimiento alternativo para el cálculo de los valores de {θ` }M
`=2
porque calcular directamente las M − 2 soluciones de la ecuación (2.139) es difı́cil.
A continuación mostramos algunos de los valores numéricos que la solución vectorial U (m, n)
y la solución exacta u(x, t) (evaluada en u(mh, nk)), toman en algunos puntos del mallado
cuando t = 1 segundo:
1
• Valor de la solución en el nodo x1 = 10 para t = 1 segundo:
   
0,03678733711496 0,03678794411714
1
U (1, 58824) =  0,03678733711496  , u , 1 =  0,03678794411714  ,
10
0,03678733711496 0,03678794411714
4
   
0,14714934845983 0,14715177646858
2
U (4, 58824) =  0,14714934845983  , u , 1 =  0,14715177646858  ,
5
0,14714934845983 0,14715177646858
6
   
0,22072402268975 0,2207276647028654
3
U (6, 58824) =  0,22072402268975  , u , 1 =  0,2207276647028654  ,
5
0,22072402268975 0,2207276647028654
9
   
0,33108603403462 0,3310914970542981
9
U (9, 58824) =  0,33108603403462  , u , 1 =  0,3310914970542981  .
10
0,33108603403462 0,3310914970542981
Para comparar numéricamente ambas soluciones vectoriales calcularemos el error global :
kU (m, n) − u(mh, nk)k∞ , y kU (m, n) − u(mh, nk)k2 .
Los errores globales para cada uno de los nodos considerados anteriormente son:

U (1, 58824) − u 1 , 1 = 0,06070021862092·10−5 ,

10
∞

U (4, 58824) − u 2 , 1 = 0,24280087448370·10−5 ,

5
∞

U (6, 58824) − u 3 , 1 = 0,36420131172832·10−5 ,

5
∞

U (9, 58824) − u 9 , 1 = 0,54630196759109·10−5 ;

10
∞
y para la norma 2

U (1, 58824) − u 1 , 1 = 0,10513586268198·10−5 ,

10
2

U (4, 58824) − u 2 , 1 = 0,42054345072791·10−5 ,

5
2

U (6, 58824) − u 3 , 1 = 0,63081517609668·10−5 ,

5
2

U (9, 58824) − u 9 , 1 = 0,94622276414261·10−5 .

10
2
m
Considerando los errores globales en todos los nodos, i.e., en xm = 10 , 1 ≤ m ≤ 9, podemos
decir que:
m
U (m, 58824) − u , 1 ≤ 0,54630196759109·10−5 , 1 ≤ m ≤ 9;

10 ∞
m
U (m, 58824) − u , 1 ≤ 0,94622276414261·10−5 , 1 ≤ m ≤ 9,

10 2
i.e., la precisión obtenida es de 4 dı́gitos.
Caso 2. M = 20
Para calcular numéricamente la solución vectorial U (m, n) de este caso, necesitaremos cono-
cer los 19 valores propios de la matriz tridiagonal Ã, pues ahora será de dimensión 19 × 19;
averiguar los valores de {θ` }19
`=2 , definidos en (2.143); y elegir un nuevo valor para el paráme-
tro r, ya que éste depende de M y de θ2 (solución de la ecuación (2.139)).

ρ2
Para este caso se tiene que − r = 4 sen θ22 = 0,05288712479686 y que el parámetro r
2
debe ser simultáneamente
r < 7,23046661988181·10−4 y r < 800 .
Tomaremos r = 7,2·10−4 . Por tanto, el paso temporal k será k = rh2 = 1,8·10−6 . Notar
que al disminuir el valor de k se necesitará un mayor número de iteraciones, n, para llegar
al mismo momento temporal t = 1.
Se muestra a continuación algunos de los valores obtenidos para este caso.
Para t = 1 los valores de U (m, n) y u(x, t) respectivamente, para algunos nodos del mallado,
son:
   
0,01839394078811 0,0183939720585721
1
U (1, 555556) =  0,01839394078811 , u , 1 =  0,0183939720585721  ,
20
0,01839394078811 0,0183939720585721
   
0,11036364472864 0,1103638323514327
3
U (6, 555556) =  0,11036364472864 , u , 1 =  0,1103638323514327  ,
10
0,11036364472864 0,1103638323514327
   
0,16554546709296 0,1655457485271491
9
U (9, 555556) =  0,16554546709296 , u , 1 =  0,1655457485271491  ,
20
0,16554546709296 0,1655457485271491
   
0,27590911182161 0,2759095808785818
3
U (15, 555556) =  0,27590911182161 , u , 1 =  0,2759095808785818  ,
4
0,27590911182161 0,2759095808785818
   
0,34948487497404 0,3494854691128702
19
U (19, 555556) =  0,34948487497404 , u , 1 =  0,3494854691128702  .
20
0,34948487497404 0,3494854691128702
Los errores globales en cada uno de estos nodos son:

U (1, 555556) − u 1 , 1 = 0,03127046500559·10−6 ,

20
∞

U (6, 555556) − u 3 , 1 = 0,18762279004048·10−6 ,

10
∞

U (9, 555556) − u 9 , 1 = 0,28143418503990·10−6 ,

20
∞

U (15, 555556) − u 3 , 1 = 0,46905697509425·10−6 ,

4
∞

U (19, 555556) − u 19 , 1 = 0,59413883518600·10−6 .

20
∞
y para la norma 2

U (1, 555556) − u 1 , 1 = 0,00541620341660·10−5 ,

20
2

U (6, 555556) − u 3 , 1 = 0,03249722050079·10−5 ,

10
2

U (9, 555556) − u 9 , 1 = 0,04874583074758·10−5 ,

20
2

U (15, 555556) − u 3 , 1 = 0,08124305125078·10−5 ,

4
2

U (19, 555556) − u 19 , 1 = 0,10290786492919·10−5 .

20
2
m
decir que:
m
U (m, 555556) − u , 1 ≤ 0,59413883518600·10−6 , 1 ≤ m ≤ 19 ;

20 ∞
m
U (m, 555556) − u , 1 ≤ 0,10290786492919·10−5 , 1 ≤ m ≤ 19 ,

20 2
i.e., la precisión obtenida para la norma infinito es de 5 dı́gitos.
Observamos numéricamente, que a medida que hacemos el mallado más fino, i.e. aumenta-
mos el valor de M , conseguimos mayor precisión.
Caso 3. M = 40
Para calcular numéricamente la solución vectorial U (m, n) de este caso, necesitaremos cono-
cer los 39 valores propios de la matriz tridiagonal Ã, pues ahora será de dimensión 39 × 39;
averiguar los valores de {θ` }39
`=2 , definidos en (2.143); y elegir un nuevo valor para el paráme-
tro r.

Para este caso se tiene que − ρr2 = 4 sen2 θ22 = 0,01292813159293 y que el parámetro r
debe ser simultáneamente
r < 8,858563511939551·10−4 y r < 3200 .
Tomaremos r = 8,8·10−4 . Por tanto, el paso temporal k será k = 5,5·10−7 .
Para t = 1 los valores de U (m, n) y u(x, t), para algunos nodos del mallado, son:
   
0,00919698258114 0,0091969860292861
1
U (1, 1818182) =  0,00919698258114 , u , 1 =  0,0091969860292861  ,
40
0,00919698258114 0,0091969860292861
   
0,09196982581135 0,0919698602928606
1
U (10, 1818182) =  0,09196982581135 , u , 1 =  0,0919698602928606  ,
4
0,09196982581135 0,0919698602928606
   
0,19313663420384 0,1931367066150072
21
U (21, 1818182) =  0,19313663420384 , u , 1 =  0,1931367066150072  ,
40
0,19313663420384 0,1931367066150072
   
0,25751551227179 0,2575156088200096
7
U (28, 1818182) =  0,25751551227179 , u , 1 =  0,2575156088200096  ,
10
0,25751551227179 0,2575156088200096
   
0,32189439033973 0,3218945110250120
7
U (35, 1818182) =  0,32189439033973 ,u , 1 =  0,3218945110250120  ,
8
0,32189439033973 0,3218945110250120
   
0,35868232066427 0,3586824551421563
39
U (39, 1818182) =  0,35868232066427 ,u , 1 =  0,3586824551421563  .
40
0,35868232066427 0,3586824551421563
Los errores globales en t = 1 son:

U (1, 1818182) − u 1 −6
,1
= 0,00344815085339·10 ,
40 ∞

U (10, 1818182) − u 1 , 1 = 0,03448150853735·10−6 ,

4
∞

U (21, 1818182) − u 21 , 1 = 0,07241116792289·10−6 ,

40
∞

U (28, 1818182) − u 7 , 1 = 0,09654822391569·10−6 ,

10
∞

U (35, 1818182) − u 7 −6
= 0,12068527982523·10 ,
,1
8 ∞

U (39, 1818182) − u 39 , 1 = 0,13447788327348·10−6 ;

40
∞
y para la norma 2

U (1, 1818182) − u 1 , 1 = 0,00597237247023·10−6 ,

40
2

U (10, 1818182) − u 1 , 1 = 0,05972372470832·10−6 ,

4
2

U (21, 1818182) − u 21 , 1 = 0,12541982187785·10−6 ,

40
2

U (28, 1818182) − u 7 −6
,1
= 0,16722642920252·10 ,
10 2

U (35, 1818182) − u 7 , 1 = 0,20903303638296·10−6 ,

8
2

U (39, 1818182) − u 39 , 1 = 0,23292252632398·10−6 .

40
2
m
decir que:
m
U (m, 555556) − u , 1 ≤ 0,13447788327348·10−6 , 1 ≤ m ≤ 39 ;

40 ∞
m
U (m, 555556) − u , 1 ≤ 0,23292252632398·10−6 , 1 ≤ m ≤ 39 ,

40 2
i.e., la precisión obtenida es de 5 dı́gitos.
Notemos que la precisión obtenida para M = 40 es mayor o igual que para los casos M = 10
y M = 20. Pero, a cambio de aumentar en exactitud debemos realizar un mayor número
de iteraciones, n, para calcular la solución en el mismo instante de tiempo t = 1. Además
también se ha observado numéricamente que el primer sumando de la solución aproximada:
M −1 n
X
2 θ` B
I − r 4 sen A− 2 sen (mθ` ) d` , (2.144)
2 M
`=2
es de magnitud despreciable al sumarse con el otro término de la solución

r n
m I + 2 B d1 , (2.145)
M
por tanto, el valor numérico de la solución aproximada construida nos lo da el término
(2.145). Pero la expresión (2.144) nos proporciona una cota para el parámetro de estabilidad
r que hace que la solución numérica construida sea una buena aproximación a la solución
exacta.
2.2. Caso general

0.4pt0.4pt 0pt0.4pt
2.2.1. Introducción
En este capı́tulo trabajaremos con sistemas de ecuaciones en derivadas parciales difusivos
fuertemente acoplados del tipo:
ut (x, t) − Auxx (x, t) − Bu(x, t) = 0, 0 < x < 1, t > 0, (2.146)

A1 u(0, t) + B1 ux (0, t) = 0, t > 0, (2.147)
A2 u(1, t) + B2 ux (1, t) = 0, t > 0, (2.148)
u(x, 0) = F (x), 0 ≤ x ≤ 1, (2.149)
donde u = (u1 , . . . , us )T y F = (f1 , . . . , fs )T son vectores s-dimensionales, elementos de Cs ,

y Ai , Bi , para i = 1, 2 son matrices complejas s × s, elementos de Cs×s verificando ciertas
propiedades a determinar.
Supondremos que

A1 B1
A= y A1 son matrices invertibles , (2.150)
A2 B2
véanse las observaciones 1.6 y 1.7.

62 2.2 Caso general
En este capı́tulo la construcción de soluciones numéricas discretas convergentes del pro-

blema (2.146)–(2.150) se realizará de la misma forma que en el apartado 2.1, es decir,
utilizando esquemas en diferencias finitas, aplicando un método de separación de variables
discreto y resolviendo explı́citamente el problema mixto (2.146)–(2.149) discretizado. Es
importante hacer notar, que el método de separación de variables discreto aquı́ propuesto,
evita la resolución de sistemas algebraicos de gran tamaño con bloques matriciales como
occurre al utilizar los métodos en diferencias estándares.
En el apartado 2.2.2 estudiaremos la discretización del problema (2.146)–(2.149) bajo la
hipótesis (2.150) y trataremos la existencia y construcción de soluciones no triviales para
el problema de contorno (2.146)–(2.149) discretizado. En el apartado 2.2.3 construiremos
soluciones discretas convergentes para el problema mixto (2.146)–(2.149) discretizado y
utilizaremos el método de las proyecciones, visto en el apartado 2.1, para extender los
resultados obtenidos a una clase más amplia de funciones de valores iniciales F (x). Final-
mente en el apartado 2.2.4 desarrollaremos dos ejemplos a modo de ilustración.
Citaremos ahora el lema 1.2 puesto que es un resultado algebraico que jugará un papel
importante en posteriores apartados de este capı́tulo y cuya demostración se encuentra en
el apartado 1.2:
Sean M y N matrices en Cs×s , entonces:
Nuc M ∩ Nuc N
h i† h i
† † †
= Im I −M M I − N I −M M N I −M M .
2.2.2. El problema de contorno discreto

Dividimos el dominio [0, 1] × [0, ∞[ en rectángulos iguales de lados ∆x = h y ∆t = k, intro-
ducimos las coordenadas de un punto tı́pico del mallado (mh, nk) y denotamos U (m, n) =
u(mh, nk). Aproximamos la derivada parcial de segundo orden utt mediante diferencias
centradas y ut mediante diferencias progresivas (véase el apartado 2.1.2) :
U (m, n + 1) − U (m, n)

ut (mh, nk) ≈
k 
, (2.151)

U (m + 1, n) − 2U (m, n) + U (m, n − 1) 
uxx (mh, nk) ≈
h2
sustituyendo (2.151) en (2.146)–(2.149) y denotando
k 1
r= , h= , (2.152)
h2 M
obtenemos la discretización del sistema de ecuaciones en derivadas parciales (2.146)–(2.149):

U (m, n + 1)
, (2.153)
rB

= rA [U (m + 1, n) + U (m − 1, n)] + I + − 2rA U (m, n) M2
1 ≤ m ≤ M − 1, n ≥ 0,
A1 U (0, n) + M B1 [U (1, n) − U (0, n)] = 0 , n≥0 (2.154)

A2 U (M, n) + M B2 [U (M, n) − U (M − 1, n)] = 0 , n≥0 (2.155)
U (m, 0) = F (mh) = f (m) , 0 ≤ m ≤ M (2.156)
Como se ha visto en el apartado 2.1.4, el esquema en diferencias (2.153) es consistente con

la ecuación (2.146) en el sentido de la definición 1.13. Mediante un método de separación
de variables discreto buscamos soluciones no triviales {U (m, n)} del problema de contorno
(2.153)–(2.155) de la forma
U (m, n) = G(n) H(m) , G(n) ∈ Cs×s , H(m) ∈ Cs . (2.157)

Sustituyendo (2.157) en (2.153) y teniendo en cuenta el apartado 2.1.3 obtenemos que
{U (m, n)} dada por (2.157) satisface (2.153) si {G(n)}, {H(m)} verifican

rB
G(n + 1) − I + 2 + ρ A G(n) = 0 , n ≥ 0 , (2.158)
M

2r + ρ
H(m + 1) − H(m) + H(m − 1) = 0 , 1 ≤ m ≤ M − 1, (2.159)
r
donde ρ es un número real. Observar que la solución de (2.158) verificando G(0) = I, viene
dada por n
rB
G(n) = I + 2 + ρ A , n ≥ 0. (2.160)
M
Si ρ verifica
−4r < ρ < 0 , (2.161)
entonces la ecuación algebraica

2r + ρ
z2 − z + 1 = 0, (2.162)
r
tiene dos soluciones diferentes z0 , z1 dadas por
2 12


2r+ρ 2r+ρ
z0 = 2r +i 1− 2r = eiθ , 



2 21 .

(2.163)
2r+ρ 2r+ρ
z1 = 2r −i 1− 2r = e−iθ , 



2r+ρ θ
ρ = −4r sen2 i2 = −1

0 < θ < π, cos θ = 2r , 2 ,
64 2.2 Caso general
Ya que la ecuación vectorial (2.159) tiene coeficientes escalares, su solución puede escribirse
de la forma
H(m) = cos(mθ) c + sen(mθ) d , c, d ∈ Cs , 1 ≤ m ≤ M − 1. (2.164)
Observación 2.4 Por (2.150) sabemos que A1 no es singular, y puesto que los bloques A2 ,
B1 , B2 pueden ser o no singulares, si en el problema de frontera discreto (2.153)–(2.155) la
matriz Ai (respectivamente Bi ) es invertible, premultiplicando la correspondiente condición
de contorno por A−1i (respectivamente por Bi−1 ) el problema queda simplificado a estudiar
los siguientes 4 casos para las condiciones de contorno:

I B1 A1 I A1 B1 A1 B1
I) , II) , III) , IV) ,
A2 B2 A2 B2 I B2 A2 I
puesto que los otros casos posibles, i.e., cuando 3 ó 4 de los bloques de A son invertibles,
son casos particulares de los cuatro arriba mencionados. Por tanto, suponer que Ai = I
para i = 1, ó 2, o Bi = I para i = 1 ó 2, no involucra pérdida de generalidad. Notemos que
en el caso I), por motivos de simplicidad, la matriz B1 representa en realidad a la nueva
matriz A−1
1 B1 6= I. Con el resto de los casos sucederı́a su correspondiente análogo.
Suponiendo que A1 = I, de (2.157) obtenemos que la condición de contorno (2.154) toma

la forma
G(n) H(0) + M B1 G(n) [H(1) − H(0)] = 0 , n ≥ 0 . (2.165)
Utilizando la ecuación (2.159) para m = 1, la ecuación (2.164) con m = 1, 2 y que cos θ =
2r+ρ
2r (al ser z0 y z1 raı́ces complejas unitarias), obtenemos H(0) = c y considerando (2.165)
para n = 0, se sigue que
[I − (1 − cos θ) M B1 ] c = −(M sen θ) B1 d . (2.166)
Premultiplicando (2.164) por [I − (1 − cos θ) M B1 ] y teniendo en cuenta (2.166) obtenemos
[I − (1 − cos θ) M B1 ] H(m)



= −M B1 cos(mθ) sen θ d + sen(mθ) [I − (1 − cos θ) M B1 ] d 

. (2.167)

θ 2m−1

= sen(mθ)I − 2M B1 sen 2 cos 2 θ d 


1≤m≤M −1 .
Por el teorema de la aplicación espectral, véase el teorema 1.6, los valores propios de la
matriz I − (1 − cos θ) M B1 son {1 − (1 − cos θ) M w ; w ∈ σ(B1 )} y la parte real de esos
valores propios son
1 − (1 − cos θ) M w1 ; w = w1 + i w2 ∈ σ(B1 ) .
Si w1 ≤ 0 entonces 1 − (1 − cos θ) M w1 6= 0. Si w1 > 0 , tomando

1
M> ,
(1 − cos θ)w1
obtenemos 1 − (1 − cos θ) M w1 < 0. Por tanto, tomando M suficientemente grande de

manera que
1
M> , (2.168)
(1 − cos θ) γ(B1 )
donde

 mı́n {Re(w) : w ∈ σ(B1 )} , si ∃w ∈ σ(B1 ), Re(w) > 0
γ(B1 ) = (2.169)
(1 − cos θ)−1 , si Re(w) ≤ 0 ∀w ∈ σ(B1 )

obtenemos que
I − (1 − cos θ) M B1 es invertible , (2.170)

y entonces para 1 ≤ m ≤ M − 1

θ 2m − 1
H(m) = sen(mθ) I − 2M sen cos θ B1 d , (2.171)
2 2
s
2 θ
solución de la ecuación (2.159) para cualquier vector d ∈ C . Notar que ρ =
es también una
−4r sen 2 para θ ∈ ]0, π[, puesto que si θ = {0, π} entonces H(m) = 0, 1 ≤ m ≤ M − 1
y por (2.157) obtendrı́amos U (m, n) = G(n) · 0 = 0, i.e., la solución trivial.
Teniendo en cuenta (2.159) para m = 1, (2.171) para m = 1, 2, que cos θ = 2r+ρ 2r y (2.165),
obtenemos
H(0) = −(M sen θ) B1 d . (2.172)
Sustituyendo (2.160), (2.171) y (2.172) en (2.165), para n > 0 se tiene
n n
rB rB
−M sen θ I + 2 + ρ A B1 − B1 I + 2 + ρ A d = 0 , n > 0. (2.173)
M M
Ya que sen θ 6= 0 para θ ∈ ]0, π[, por (2.163) tenemos

r −1
w= 2
= θ
6= 0 , (2.174)
M ρ 4M 2 sen2 2
y (2.173) puede reescribirse de la siguiente forma
[(I + ρ (A + w B))n B1 − B1 (I + ρ (A + w B))n ] d = 0 , d ∈ Cs , n > 0 . (2.175)
Considerando (2.159) para m = M − 1, es fácil probar que

θ 2M − 1
H(M ) = sen(M θ) I − 2M sin cos θ B1 d , d ∈ Cs . (2.176)
2 2
Imponiendo a U (m, n), dada por (2.157), la condición de contorno (2.155) para n ≥ 0 y
haciendo uso de (2.160), (2.171) y (2.176) obtenemos
{A2 (I + ρ(A + wB))n sen(M θ)

66 2.2 Caso general

θ 2M − 1
−2M sen cos θ A2 (I + ρ(A + wB))n B1
2 2

θ 2M − 1
+2M sen cos θ B2 (I + ρ(A + wB))n
2 2

2 2 θ n
+ 4M sen sen((M − 1)θ)B2 (I + ρ(A + wB)) B1 d = 0 , n ≥ 0. (2.177)
2
Sustituyendo (2.175) en (2.177) para n > 0 y utilizando (2.177) para n = 0, se sigue que
para n ≥ 0

θ 2M − 1
A2 sen(M θ) − 2M sen cos θ A2 B1
2 2

θ 2M − 1
+2M sin cos θ B2
2 2

θ
+4M 2 sen2 sen((M − 1)θ)B2 B1 (I + ρ(A + wB))n d = 0. (2.178)
2
Sea p el grado del polinomio minimal de A + w B, entonces por el teorema de Cayley-
Hamilton, véase el teorema 1.7, para n ≥ p las potencias (A + w B)n se expresan en
términos de I, A + w B, (A + w B)2 , . . . , (A + w B)p−1 . Ya que w 6= 0, las condiciones
(2.175) y (2.178) se cumplen si:
[(A + wB)n B1 − B1 (A + wB)n ] d = 0 , 0 < n < p. (2.179)

y

θ 2M − 1
2M sen cos θ (B2 − A2 B1 ) + A2 sen(M θ)
2 2

2 2 θ
+4M sen sen((M − 1)θ)B2 B1 (A + wB)n d = 0 , (2.180)
2
0 ≤ n < p.
respectivamente. Para garantizar que {U (m, n)} sea una solución no trivial, el vector d que
aparece en (2.180) debe ser no nulo. Por (2.180), existirá un vector no nulo d verificando
la condición (2.180) si la matriz

θ 2M − 1
L(θ) = 2M sen cos θ (B2 − A2 B1 ) + A2 sen(M θ)
2 2

2 2 θ
+4M sen sen((M − 1)θ) B2 B1 es singular, (2.181)
2
0 < θ < π.
Observemos que sumando y restando el término A2 sen((M − 1)θ) la matriz L(θ) puede
reescribirse de la forma:

θ 2M − 1 A2
L(θ) = 2M sen cos θ (B2 − A2 B1 ) +
2 2 M

θ
+ sen((M − 1)θ) A2 + 4M 2 sen2 B2 B1 . (2.182)
2
Por las propiedades del complemento de Schur de una matriz, véasen las observaciones 1.6
y 1.7, junto con la hipótesis (2.150) y A1 = I, se sigue que
B2 − A2 B1 es invertible. (2.183)
Por (2.183) y el lema de Banach, véase el lema 1.3, se sigue que
A2
(B2 − A2 B1 ) + es invertible si M > kA2 k (B2 − A2 B1 )−1 . (2.184)

M
Por tanto si M satisface la condición (2.184), los valores de θ ∈ ]0, π[ para los cuales la
matriz L(θ) definida por (2.182) es singular deberán verificar que sen((M − 1)θ) 6= 0.
Consecuentemente L(θ) es singular si y sólo si

2 2 θ
A2 + 4M sen B2 B1 +
2
2M sen 2θ cos 2M2−1 θ

A2
+ (B2 − A2 B1 ) + es singular, (2.185)
sen ((M − 1)θ) M
o equivalentemente

sen(M θ) −1 2 θ
2
(B2 − A2 B1 ) A2 + 4M sen (B2 − A2 B1 )−1 B2 B1 +
sen((M − 1)θ) 2
2M sen 2θ cos 2M2−1 θ

I , es singular, 0 < θ < π. (2.186)
sen ((M − 1)θ)
Para simplificar la expresión (2.186) introducimos las matrices
b2 = (A2 B1 − B2 )−1 A2 ,
A b2 = (A2 B1 − B2 )−1 B2 = A
B b2 B1 − I . (2.187)
Utilizando las matrices A b2 , B

b2 definidas en (2.187) y el teorema de la aplicación espectral
la condición (2.186) significa que

sen(M θ)
M −1 es un valor propio de la matriz
sen ((M − 1)θ) 

. (2.188)

sen(M θ) b2 + 4M 2 sen2 θ

b2 B 2 − B1 ,

A 2 A 1 0<θ<π
sen ((M − 1)θ)
Supongamos que

Existen α ∈ σ A b2 ∩ R ; β ∈ σ(B1 ) ∩ R y v ∈ Cs ∼ {0}
. (2.189)
tales que A b2 − α I v = (B1 − β I) v = 0
Por (2.189) se sigue que

68 2.2 Caso general

sen(M θ) b2 + 4M 2 sen2 θ b2 B 2 − B1 v =
A A 1
sen ((M − 1)θ) 2

sen(M θ) 2 2 θ 2

= α + 4M sen αβ − β v , 0 < θ < π,
sen ((M − 1)θ) 2
o equivalentemente

v es un vector propio de la matriz


sen(M θ) b2 + 4M 2 sen2 θ b2 B 2 − B1

A A 1

sen ((M − 1)θ) 2 

. (2.190)

asociado al valor propio real 



sen(M θ) θ 
α + 4M 2 sen2 αβ 2 − β

sen ((M − 1)θ) 2
Tomando M suficientemente grande de manera que
M > α,
la condición (2.188) y (2.190) hacen posible encontrar soluciones de la ecuación escalar
αβ 2 − β

sen(M θ) M 2 2 θ
= + 4M sen , 0 < θ < π,
sen((M − 1)θ) M −α 2 M −α
o equivalentemente
cot ((M − 1) θ)

M 1 2
θ
= − cot θ + + 2M αβ − β tan , 0 < θ < π. (2.191)
M − α sen θ 2
i h
Para cada entero ` con 1 ≤ ` ≤ M − 1, en el intervalo J` = (`−1)π `π
M −1 M −1 se satisface
,
lı́mθ→ (`−1)π + cot((M − 1)θ) = +∞ ;


M −1


lı́mθ→ `π − cot((M − 1)θ) = −∞ ;  (2.192)

M −1 

cot((M − 1)θ) decreciente en J` ,
donde la última condición se cumple porque
d M −1
(cot((M − 1)θ)) = − 2
< 0.
dθ sen ((M − 1)θ)
Además la función eM (θ) describiendo la parte derecha de la igualdad (2.191) es continua

y creciente en ]0, π[ si se satisface alguna de las siguientes condiciones

β = 0,
ó 

αβ = 1 ,



α<1 y ó .

(2.193)
β > 0 y αβ > 1 , 


ó 

β < 0 y αβ < 1 .
Entonces por (2.192)–(2.193) existe una única solución θ` de (2.191) en el intervalo J` ,

verificando

((M − 1) θ` )
cot
M 1 2
θ` 
= − cot θ` + + 2M αβ − β tan 
.
M − α sen θ` 2  (2.194)

1 ≤ ` ≤ M − 1, θ` ∈ J`
Por tanto la condición (2.180) puede escribirse de la forma
S(α, β, θ` ) (A + w` B)n d` = 0 , (2.195)

0 ≤ n ≤ p` − 1 , 1 ≤ ` ≤ M − 1 ,
donde

sen(M θ` ) 2 2 θ`

b2 B 2 − B1 +
S(α, β, θ` ) = A2 + 4M sen
b A 1
sen ((M − 1)θ` ) 2

sen(M θ` ) 2 2 θ` 2

− α + 4M sen αβ − β I ,
sen ((M − 1)θ` ) 2
(2.196)
p` es el grado del polinomio minimal de la matriz A + w` B, siendo θ` la solución de (2.194)

y
−1
w` = , 1 ≤ ` ≤ M − 1. (2.197)
θ`
4M 2 sen2
2
Además en vista de (2.194) la ecuación (2.179) toma la forma
[(A + w` B)n B1 − B1 (A + w` B)n ] d` = 0 , 0 < n < p` . (2.198)

A fin de determinar los vectores d` ∈ Cs ∼ {0} que verifican las condiciones (2.195) y
(2.198) introducimos la matriz por bloques definida por
70 2.2 Caso general
 
B1 (A + w` B) − (A + w` B)B1
 B1 (A + w` B)2 − (A + w` B)2 B1 
..
 
 
 . 
 B1 (A + w` B)p` −1 − (A + w` B)p` −1 B1
 

 
T (α, β, θ` ) = 
 S(α, β, θ` ) .
 (2.199)

 S(α, β, θ` )(A + w` B) 


 S(α, β, θ` )(A + w` B)2 

 .. 
 . 
S(α, β, θ` )(A + w` B)p` −1
Entonces utilizando la definición de T (α, β, θ` ) los vectors d` deberán satisfacer
T (α, β, θ` )d` = 0 , 1 ≤ ` ≤ M − 1, d` ∈ Cs ∼ {0} . (2.200)
La ecuación (2.200) tendrá soluciones no nulas d` ∈ Cs si sólo si
rango (T (α, β, θ` )) < s .

Pero no vamos a utilizarla porque es muy indirecta. Ası́ que buscaremos condiciones su-
−1
ficientes más fáciles de comprobar que nos garanticen que {d` }M
`=1 cumplen la ecuación
(2.200).
−1
Imponiendo que los vectores {d` }M
`=1 satisfagan las condiciones:

(B1 − βI) d` = A b2 − αI d` = 0 , d` ∈ Cs ∼ {0} , 1 ≤ ` ≤ M − 1 , (2.201)
y
{(A + w` B)n d` ; 1 ≤ n ≤ p` − 1}

⊂ Nuc A2 − αI ∩ Nuc (B1 − βI) , 1 ≤ ` ≤ M − 1 , (2.202)
b
−1
garantizamos que los {d` }M`=1 verifiquen (2.195) y (2.198), es decir, garantizamos que
−1
{d` }M
`=1 sean soluciones no nulas de la ecuación (2.200) para cada θ` ∈ J` , 1 ≤ ` ≤ M − 1,
de la ecuación (2.194).
En efecto,
(i) [B1 (A + w` B)n − (A + w` B)n B1 ] d`

= B1 (A + w` B)n d` − (A + w` B)n B1 d`
= B1 (A + w` B)n d` − β (A + w` B)n d`
+β (A + w` B)n d` − (A + w` B)n B1 d`
= (B1 − βI)(A + w` B)n d` − (A + w` B)n (B1 − β I) d` = 0 ,
ya que de (2.201) y (2.202) tenemos que
(B1 − β I) d` = 0 y (A + w` B)n d` ∈ Nuc(B1 − β I) .

(ii) S(α, β, θ` )(A + w` B)n d`

sen(M θ` )
b2 − αI (A + w` B)n d`
= A
sen ((M − 1)θ` )

2 2 θ`
+4M sen ·
2
h i
Ab2 B 2 − B1 − (αβ 2 − β)I (A + w` B)n d`
1
sen(M θ` )
b2 − αI (A + w` B)n d`
= A
sen ((M − 1)θ` )

θ`
+4M 2 sen2 ·
2
h i
Ab2 B 2 − Ab2 β 2 + Ab2 β 2 − B1 − (αβ 2 − β)I (A + w` B)n d`
1
sen(M θ` )
b2 − αI (A + w` B)n d`
= A
sen ((M − 1)θ` )
h
2 2 θ` b2 (B1 + βI)(B1 − βI)(A + w` B)n d`
+4M sen A
2
i
+β 2 A b2 − αI (A + w` B)n d` − (B1 − βI)(A + w` B)n d` ,
0 ≤ n ≤ p` − 1 .
Para n = 0 obtenemos:
sen(M θ` )
b2 − αI d`
S(α, β, θ` ) d` = A
sen ((M − 1)θ` )

θ`
+4M 2 sen2 ·
2
h i
b2 (B1 + βI)(B1 − βI) + β 2 A
A b2 − αI − (B1 − βI) d` = 0 ,
ya que se cumple la condición (2.201). Para 1 ≤ n ≤ p` − 1 obtenemos que:

S(α, β, θ` )(A + w` B)n d` = 0 ,
ya que se cumple la condición (2.202). Reemplazando θ por θ` en (2.160) y (2.171); utili-
zando (2.157); y que d` son vectores propios de la matriz B1 , es decir (B1 − β I)d` = 0, se
sigue que
n
θ` B
U` (m, n) = I − r 4 sin2 A− 2 ·
2 M

θ` 2m − 1
· sen (mθ` ) − 2M β sen cos θ` d` , (2.203)
2 2
para 1 ≤ m ≤ M − 1, n ≥ 0, define soluciones no triviales del problema de contorno
(2.153)–(2.155).
72 2.2 Caso general
Resumiendo se ha demostrado el siguiente resultado.
TEOREMA 2.5 Consideremos el problema de contorno (2.153)–(2.155) bajo la hipóte-

sis (2.150) con A1 = I, sea M > 0 un entero positivo suficientemente grande de ma-
nera que se cumplen las condiciones (2.168) y (2.184), y consideremos la matriz A
b2 =
−1
(A2 B1 − B2 ) A2 .
(i) Asumimos
i las condiciones
h (2.189) y (2.193). Entonces existen soluciones θ` de (2.194),
(`−1)π `π
θ` ∈ M −1 , M −1 = J` , 1 ≤ ` ≤ M − 1, para las cuales la matriz L(θ` ) definida por
(2.182) es singular.
(ii) Bajo las hipótesis vistas en (i), consideremos d` vectores in Cs verificando (2.201) y
(2.202) para 1 ≤ ` ≤ M − 1, entonces {U` (m, n)} dada por (2.203) define soluciones
no triviales del problema de contorno (2.153)–(2.155).
Observación 2.5 El caso donde además de la invertibilidad de A tenemos B1 = I se trata

de forma análoga al caso estudiado, teniendo en cuenta las propiedades del complemento
de Schur, véase la observación 1.6. Por otro lado, considerando el cambio
m → M − m, 0≤m≤M,
los casos donde A2 = I o B2 = I se transforman en los casos anteriores.
2.2.3. El problema mixto discreto

En este apartado estudiaremos la construcción de soluciones exactas del problema mixto
discreto (2.153)–(2.156). Asumimos la notación y las hipótesis del teorema 2.5-(i) y (ii).
Por superposición de las M − 1 soluciones obtenidas para el problema de contorno (2.153)–
(2.155) obtenemos

U (m, n)

M −1
n 
P 2 θ` B 
= I − r 4 sen A− 2 · 
`=1 2 M 


. (2.204)
βM ρ` 
· 1− sen(mθ` ) − βM cos(mθ` ) sen(θ` ) d` , 
2r 



θ `

ρ` = −4r sen2 , 1≤`≤M −1
2
Imponiendo a {U (m, n)}, dada por (2.204), la condición (2.156) obtenemos que los vectores
d` que aparecen en (2.204) deben satisfacer
M −1
X βM ρ`
f (m) = 1+ sen(mθ` ) − βM cos(mθ` ) sen(θ` ) d` . (2.205)
2r
`=1
M −1
Bajo la hipótesis (2.150), si {H` (m)}`=1 son las soluciones del problema de Sturm-Liouville
discreto:
ρ 
−H(m + 1) + 2H(m) − H(m − 1) = − H(m) 
r







βM 

H(0) = H(1) , (2.206)
βM − 1 




M (αβ − 1) 

H(M ) = H(M − 1) 

α + M (αβ − 1)

entonces por el apartado 2.2.2 se obtiene que:

βM ρ`
H` (m) = 1+ sen(mθ` ) − βM cos(mθ` ) sen(θ` ) d` ,
2r
1≤m≤M −1
i h
M −1
donde {θ` }`=1 ∈ (`−1)π `π
M −1 , M −1 son las soluciones de la ecuación (2.194) y los vectores
−1
{d` }M
`=1 verifican las condiciones (2.201) y (2.202).
Como los coeficientes del problema (2.206) son escalares, consideramos el siguiente problema
de Sturm-Liouville discreto escalar:
ρ 
−h(m + 1) + 2h(m) − h(m − 1) = − h(m) 
r







βM 

h(0) = h(1) , 1 ≤ m ≤ M − 1. (2.207)
βM − 1 




M (αβ − 1) 

h(M ) = h(M − 1) 

α + M (αβ − 1)

tiene exactamente M − 1 valores propios dados

Por el apartado 1.2, el problema (2.207)
−ρ` M −1
por r `=1 , donde ρ` = −4r sen2 θ2` y θ` verifica (2.194). Para cada valor propio −ρ r
`
existe una función propia

βM ρ`
{h` (m)} = 1+ sen(mθ` ) − βM cos(mθ` ) sen(θ` ) , (2.208)
2r
y estas funciones propias son ortogonales en N (1, M − 1) con respecto a la función peso
w(m) = 1, para 1 ≤ m ≤ M − 1.
Puesto que f (m) y d` son vectores en Cs podemos escribir (2.205) de forma escalar:
M −1
X βM ρ`
fq (m) = 1+ sen(mθ` ) − βM cos(mθ` ) sen(θ` ) d`,q , (2.209)
2r
`=1
donde fq (m) y d`,q denotan la q-ésima componente de los vectores f (m) y d` respectiva-
mente. Por la ortogonalidad de las funciones propias {h` (m)} que aparecen en (2.205) y la
teorı́a de las series de Fourier discretas, véase el apartado 1.2, se sigue que
74 2.2 Caso general
M −1 n o
P βM ρ`
1+ 2r sen(νθ` ) − βM cos(νθ` ) sen(θ` ) fq (ν)
ν=1
d`,q = M −1 n o2 ,
P βM ρ` (2.210)
1+ 2r sen(νθ` ) − βM cos(νθ` ) sen(θ` )
ν=1
1 ≤ ` ≤ M − 1, 1 ≤ q ≤ s,
o equivalentemente en forma vectorial
M −1 n o
P βM ρ`
1+ 2r sen(νθ` ) − βM cos(νθ` ) sen(θ` ) f (ν)
ν=1
d` = M −1 n o2 ,
P βM ρ` (2.211)
1+ 2r sen(νθ` ) − βM cos(νθ` ) sen(θ` )
ν=1
1 ≤ ` ≤ M − 1.
La expresión (2.211) para vectores d` debe ser compatible con las condiciones (2.201) y
(2.202). Esto significa que {f (m)} debe verificar

(B1 − β I) f (m) = A b2 − α I f (m) = 0 , 1 ≤ m ≤ M − 1 , (2.212)
y si w` está definido en (2.197),

{(A + w` B)n f (m) , 1 ≤ n ≤ p` − 1} ⊂ Nuc Ab2 − α I ∩ Nuc (B1 − β I) , (2.213)
para 1 ≤ m ≤ M − 1, 1 ≤ ` ≤ M − 1.
−1
Si {f (m)}Mm=1 verifica (2.212) y (2.213) entonces {U (m, n)} definida por (2.204), donde
d` están dados en (2.211), es una solución del problema (2.153)–(2.156). Observar que las
condiciones (2.212) y (2.213) se satisfacen si

f (m) ∈ Nuc A b2 − α I ∩ Nuc (B1 − β I) , 1 ≤ m ≤ M − 1 , (2.214)
y
#
Nuc Ab2 − α I ∩ Nuc (B1 − β I) es un subespacio invariante
, (2.215)
por la matriz A + w` B , 1 ≤ ` ≤ M − 1.

Utilizando el lema 1.2, las condiciones (2.214) y (2.215) se reescriben de la forma
f (m) ∈ Im L(α, β) , 1 ≤ m ≤ M − 1, (2.216)


I − L(α, β)L(α, β)† (A + w` B) L(α, β) = 0 , 1 ≤ ` ≤ M − 1, (2.217)
donde
h i† h i
L(α, β) = I − Pα† Pα I − Qβ I − Pα† Pα Qβ I − Pα† Pα 

 (2.218)
b2 − α I ,
Pα = A Qβ = B1 − β I ,
Notemos que la condición (2.217) significa que Im L(α, β) es un subespacio invariante por
la matriz A + w` B, para 1 ≤ ` ≤ M − 1, véase la observación 1.3. La solución {U (m, n)}
del problema mixto discreto (2.153)–(2.156), definida por (2.204), (2.211), es estable, i.e.
permanece acotada cuando n → ∞ si {f (m)} está acotado y las matrices

2 θ` B
I − r 4 A sen − 2 , 1 ≤ ` ≤ M − 1,
2 M
son convergentes, véase el lema 1.1. Por el teorema 2.2 del apartado 2.1 esto ocurre si
A + AH

x>0 para todo x ∈ σ , (2.219)
2
B + BH

y≤0 para todo y ∈ σ , (2.220)
2
H H H H
e1 = A+A , B
y siiA e1 = B+B , A e2 = A−A , B e2 = B−B y θ1 es la única solución de (2.194)
h2 2 2i 2i
π
en 0, M −1 , r debe verificar
" 2 #
2 θ1
e1 − λmax B
M 2M sen λmı́n A e1
2
r<h i2 h i2 . (2.221)
4M 2 λmax A e1 + ρ B e1 + 4M 2 λmax Ae2 + ρ B e2
Resumiendo se ha demostrado el siguiente resultado.
TEOREMA 2.6 Consideremos el problema mixto (2.153)-(2.156) con las hipótesis (2.189)
y (2.150) con A1 = I. Sea A b2 = (A2 B1 − B2 )−1 A2 y M > 0 un número entero sufi-
cientemente grande de forma que verifique las condiciones (2.168), (2.184) y (2.193). Sea
θ` la solución de (2.194) y w` definido por (2.197) para 1 ≤ ` ≤ M − 1. Supongamos
que {f (m)} satisface la condición (2.216) y las matrices A + w` B verifican (2.217), don-
de L(α, β) está definida por (2.218). Entonces {U (m, n)} definida por (2.204), donde d`
vienen dados por (2.211), es una solución del problema mixto discreto (2.153)–(2.156).
Además, {U (m, n)} es estable si: las matrices A, B satisfacen las condiciones (2.219)–
(2.220), {f (m)} está acotado y r es suficientemente pequeño de manera que se cumple
(2.221), entonces .
76 2.2 Caso general
Ahora estudiaremos condiciones más generales que las consideradas en el teorema 2.6.
Supongamos que

Λ = {α1 , . . . , αt } ⊂ R ∩ σ A
b2 , (2.222)
Ω = {β1 , . . . , βq } ⊂ R ∩ σ(B1 ) . (2.223)

Por el lema 1.2 la condición
L (αi , βj ) 6= 0 , 1 ≤ i ≤ t, 1 ≤ j ≤ q, (2.224)
es equivalente a

Nuc A b2 − αi I ∩ Nuc (B1 − βj I) 6= ∅ , 1 ≤ i ≤ t, 1 ≤ j ≤ q. (2.225)
Consideramos el conjunto F ⊂ Λ × Ω definido por
 
 i` , βj` ) ∈
 (α Λ × Ω , cumpliendo alguna condición de (2.193) 

F= A2 − αi` I v` = (B1 − βi` I) v` = 0 , v` ∈ Cs ∼ {0} ,
b (2.226)
 
L (αi` , βj` ) 6= 0
 
y la matriz por bloques
L = L (αi1 , βj1 ) , L (αi2 , βj2 ) , . . . , L αip , βjp ∈ Cs×ps ,

(2.227)
y supongamos que f (m) ∈ Im L para 0 ≤ m ≤ M , o equivalentemente

I − LL† f (m) = 0 , 0 ≤ m ≤ M , (2.228)
porque Im L = Nuc I − LL† . Por el lema 1.2 obtenemos

Sδ = Im L (αiδ , βjδ ) = Nuc Ab2 − αi I ∩ Nuc (B1 − βj I) ,
δ δ
(2.229)
1≤δ≤p
y por (2.227), (2.229), el subespacio Im L es suma directa de los subespacios Sδ ,
Im L = S1 ⊕ S2 ⊕ · · · ⊕ Sp . (2.230)
n oM
Sea fbδ (m) la proyección de {f (m)}M
m=0 sobre el subespacio Sδ , definida por
m=0
fbδ (m) = [0, . . . , 0, L (αiδ , βjδ ) , 0, . . . , 0] L† f (m) , (2.231)

1 ≤ δ ≤ p, 0 ≤ m ≤ M .
Ya que fbδ (m) está en Sδ , por (2.228) se sigue que

p
X
fbδ (m) = LL† f (m) = f (m) , 0≤m≤M. (2.232)
δ=1
(δ)
Supongamos que Im L (αiδ , βjδ ) es un subespacio invariante por la matriz A + w` B, i.e.
h i 
(δ)
I − L (αiδ , βjδ ) L (αiδ , βjδ )† A + w` B L (αiδ , βjδ ) = 0 ,


−1

(δ) , (2.233)
w` = (δ)
! , 1 ≤ ` ≤ M − 1, 
θ `

4M 2 sen2

2
(δ)
véase la observación 1.3, donde θ` es la solución de (2.194) asociada al par (αiδ , βjδ ) en
J` .
Consideramos el problema (Pδ ) definido por (2.153)–(2.155) junto con la condición inicial
U (m, 0) = fbδ (m) , 0 ≤ m ≤ M , 1 ≤ δ ≤ p, (2.234)
y observemos que la solución {Uδ (m, n)} del problema (Pδ ) está definida por (2.204) donde
(δ)
d` viene dados por
( (δ)
! )
M −1 βjδ M ρ`
P (δ) (δ) (δ)
1− sen νθ` − βjδ M cos νθ` sen θ` fbδ (ν)
(δ) ν=1 2r
d` = ,
2 (2.235)
( (δ)
! )
M −1 βjδ M ρ`
P (δ) (δ) (δ)
1− sen νθ` − βjδ M cos νθ` sen θ`
ν=1 2r
para 1 ≤ ` ≤ M − 1 , 1 ≤ δ ≤ p , 1 ≤ j ≤ q, donde
M −1
" (δ)
! !#n
X
2 θ` B
Uδ (m, n) = I − r 4A sen − 2 ·
2 M
δ=1
" (δ)
! #
βj M ρ` (δ) (δ) (δ) (δ)
· 1− δ sen(mθ` ) − βjδ M cos(mθ` ) sen(θ` ) d` . (2.236)
2r
Por linealidad y (2.232), (2.236) se sigue que
p
X
U (m, n) = Uδ (m, n) , 1 ≤ m ≤ M − 1, n ≥ 0, (2.237)
δ=1
es una solución del problema mixto discreto (2.153)–(2.156). Además (2.237) es una solución
estable si se cumplen las condiciones (2.219)–(2.220), {f (m)} está acotado y el parámetro
r verifica
78 2.2 Caso general
  !!2  
(δ)
 θ1 
M 2  2M sen
 e1 − λmax B 

 λmı́n A e1  



 2 


r < mı́n h i2 h i2 . (2.238)
1≤δ≤p 
4M 2 λ A1 + ρ B1 2
+ 4M λmax A2 + ρ B

max
 e e e e2 


 


 

Resumiendo el siguiente resultado es una consecuencia del teorema 2.6.
TEOREMA 2.7 Consideremos el problema (2.153)–(2.156) bajo la hipótesis (2.150) con

A1 = I, supongamos que se satisfacen las condiciones (2.222) y (2.223) y sea M un núme-
ro entero verificando (2.168) y (2.184). Sean F y L definidas por (2.226) y (2.227) res-
pectivamente, supongamos que {f (m)} está acotado, que las condiciones (2.219)–(2.220)
se satisfacen y que r es suficientemente pequeño de manera que (2.238) se cumple. Sea
n oM
(δ)
fbδ (m) definida por (2.231), w definida por (2.233) y supongamos que se verifica
`
m=0
la condición (2.233). Si {Uδ (m, n)} viene dada por (2.236) entonces {U (m, n)} definida
por (2.237) es una solución estable del problema mixto discreto (2.153)–(2.156) .
2.2.4. Ejemplos
EJEMPLO 2.4 Consideramos el problema (2.146)–(2.149) con los datos:
   1 
1 0 1 −2 0 1
   
1 1
   
A= 0 2 , B =  0 − 4
0  0 , A1 = I ,

   
1
−1 0 2 −1 0 −6
1
− 12 4
     
2 0 1 1 0 1 0 3
     
1
     
 1 −1  0 −1 1
B1 =  1 , A2 = 
 , B2 =  −1 0 ,
  3 
     
1 1 1
0 0 3 1 0 1 2 0 3
f (m) = F (mh) = (f1 (m), f2 (m), f3 (m))T , y h = 1

M , 1 ≤ m ≤ M − 1.
La hipótesis (2.150) se satisface, y además tenemos A b2 = (A2 B1 − B2 )−1 A2 = A2 con

1 1
σ A2 = {−1, 0, 2} ,
b σ (B1 ) = −1, , .
3 2
Sean α1 = −1, α2 = 0, α3 = 2, β1 = −1 y β2 = 21 . Los pares de valores propios que verifican

al menos una de las condiciones dadas en (2.193) son
(α1 , β1 ) = (−1, −1) ,

(α2 , β1 ) = (0, −1) ,

(α3 , β1 ) = (2, −1) ,

1
(α3 , β2 ) = 2, .
2
Ahora debemos seleccionar aquellos pares de valores propios que tengan un vector propio
común asociado a ellos.
   
0 0
b2 − α1 I
A v1 = 0 ⇐⇒ v1 =  v1,2  = v1,2  1  ;
0 0
    
v2,1 1 0
b2 − α2 I
A v2 = 0 ⇐⇒ v2 = −v2,1 = v2,1 −1 + v2,3 0  ;
    
v2,3 0 1
     
v3,1 1 0
b2 − α3 I
A v3 = 0 ⇐⇒ v3 =  13 v3,1  = v3,1  13  + v3,3  0  ;
v3,3 0 1
   
0 0
b1 − β1 I
B w1 = 0 ⇐⇒ w1 = w1,2 = w1,2 1  ;
  
0 0
   
w2,1 1
b1 − β2 I
B w2 = 0 ⇐⇒ w2 =  23 w1,2  = w2,1  32  .
0 0
 
0
Por tanto ∃ v =  1 , tal que
0

b2 − α1 I
A v = (B1 − β1 I) v = 0 ,
y sólo podremos trabajar con el par de valores propios (α1 , β1 ) = (−1, −1). Ahora debemos
asegurarnos de que las matrices A + w` B, 1 ≤ ` ≤ M − 1, satisfacen la condición (2.217),
i.e.:
h i
I − L1 (α1 , β1 ) L1 (α1 , β1 )† (A + w` B) L1 (α1 , β1 ) = 0 , 1 ≤ ` ≤ M − 1 ,
donde
h i† h i
L1 (α1 , β1 ) = I − Pα†1 Pα1 I − Qβ1 I − Pα†1 Pα1 Qβ1 I − Pα†1 Pα1 ,
b2 − α1 I , Qβ = B1 − β1 I .
Pα 1 = A 1
80 2.2 Caso general
Para nuestro caso obtenemos:
4
− 72 − 17 3
     
2 0 1 7 2 0 1
     
 , Pα† = 
     
Pα1  0 0 1
= 0 0 0  , Qβ1 = 
 1 0 0 
;

 1 
     
3 5 3 4
1 0 2 − 14 14 7 0 0 3
 
0 0 0
  h i†
Pα†1 Pα1  , Qβ I − Pα† Pα = Qβ I − Pα† Pα
 
I−  0 1 0
=  1 1 1 1 1 1 =0 .
 
0 0 0
Por tanto tenemos:

 
0 0 0
6 0
L (α1 , β1 ) =  0 1 0  = 

0 0 0 

. (2.239)
  
1 0 0 
†

I − L (α1 , β1 ) L (α1 , β1 ) =  0 0 0  
0 0 1
Además  w` 
1− 2 0 1 + w`
 
1 w`
 
A + w` B = 
 0 2 − 4 0 .
 (2.240)
 
w`
−1 − w` 0 2− 6
Por (2.239) y (2.240) se sigue que
h i
I − L(−1, −1)L(−1, −1)† (A + w` B) L(−1, −1) = 0 , (2.241)
1 ≤ ` ≤ M − 1,
La condición (2.228) se satisface para cualquier función vectorial {f (m)} de la forma
f (m) = (0, f2 (m), 0))T .
Observamos que
 
1 0 0
A + AH A + AH
 
 1
 1
 0
= 0 , σ = , 1, 2 ;

2 2 2 2
 
0 0 2
− 21
 
0 0
B + BH B + BH
 
1 1 1
− 14
 
= 0 0 , σ
 = − ,− ,− ,
2 
  2 2 4 6
0 0 − 16
y por tanto las condiciones de estabilidad (2.219)–(2.220) se satisfacen. Tomando r sufi-
cientemente pequeño de forma que verifique la condición (2.238), y M verificando (2.168)
y (2.184), por el teorema 2.7 se obtiene que la sucesión vectorial

U (m, n)

M −1
n 
P 2 θ` B 
= I − r 4 sen A− 2 · 
`=1 2 M 


,
β1 M ρ` 
· 1+ sen(mθ` ) − β1 M cos(mθ` ) sen(θ` ) d` , 
2r 



θ `

ρ` = −4r sen2 , 1≤`≤M −1
2
es una solución numérica exacta del problema (2.153)–(2.156) donde

M −1 n o
P M β 1 ρ`
1+ 2r sen(νθ` ) − M β1 cos(νθ` ) sen(θ` ) f (ν)
ν=1
d` = M −1 n o2 .
P M β1 ρ `
1+ 2r sen(νθ` ) − M β1 cos(νθ` ) sen(θ` )
ν=1
1≤`≤M −1
−1
y {θ` }M
`=1 son las soluciones de la ecuación (2.194).
82 2.2 Caso general
EJEMPLO 2.5 Consideremos el problema (2.146)–(2.149) con los datos:

   
1 −1 0 −3 2 0
   
   
A=  0 2 0 , B =  0 −8
  0 
, A1 = I ,
   
0 2 1 0 5 −3
3
   
3 1 0 
2 1 0

   5 0
2 
     
 1     
B1 =  0 − 2 0 , A2 =  0 2 0 
, B2 =  0 −2 0 ,
    
     
1
     
 1 
0 −1
 9 1 
0 2 − 2 0 − −
2 4 2
f (m) = F (mh) = (f1 (m), f2 (m), f3 (m))T , y h = 1
M , 1 ≤ m ≤ M − 1.
b2 = (A2 B1 − B2 )−1 A2 = A2 con
La hipótesis (2.150) se satisface, y A

1
σ A2 = {−1, 2} ,
b σ (B1 ) = − , 3 .
2
Sean α1 = −1, α2 = 2, β1 = − 21 , β2 = 3 y observemos que ambos pares (α1 , β1 ), (α2 , β2 )
satisfacen alguna de las condiciones dadas en (2.193) y además
 
0
v= 0 , b2 − α1 I v = (B1 − β1 I) v = 0 ,
A
1
 
1
w= 0 , b2 − α2 I w = (B1 − β2 I) w = 0 .
A
0
Para el par (α1 , β1 ) = (−1, −1/2) la matriz L(α1 , β1 ) definida por (2.218) toma el valor

 
0 0 0
L (−1, −1/2) =  0 0 0  6= 0 

0 0 1 

. (2.242)
  
1 0 0 
†

I − L (−1, −1/2) L (−1, −1/2) = 0 1 0
  
0 0 0
(1) M −1
n o
Sean θ` las soluciones de (2.194) correspondientes al par (−1, −1/2) y sea
`=1
(1) −1
w` = (1)
!, 1 ≤ ` ≤ M − 1.
θ`
4M 2 sin2
2
Por tanto,
(1) (1)
 
1 − 3w` −1 + 2w` 0
 
 
(1) (1)
A + w` B =  0 2 − 8w` 0 . (2.243)
 
 
 
(1) (1)
0 2 + 5w` 1 − 3w`
Por (2.242) y (2.243) se sigue que
h i
† (1)
I − L (−1, −1/2) L (−1, −1/2) A+ w` B L (−1, −1/2) = 0 , (2.244)
1 ≤ ` ≤ M − 1.
Consideremos ahora el par (α2 , β2 ) = (2, 3). Haciendo cálculos obtenemos
   
1 0 0 0 0 0
6 0
L(2, 3) =  0 0 0  = y I − L(2, 3)L(2, 3)† =  0 1 0  .
0 0 0 0 0 1
(2) M −1
n o
Sean θ` las soluciones de (2.194) correspondientes al par (2, 3) y sea
`=1
(2) −1
w` = (2)
!, 1 ≤ ` ≤ M − 1.
θ`
4M 2 sin2
2
Notemos que
(2) (2)
 
1 − 3w` −1 + 2w` 0
 
 
(2) (2)
A + w` B =  0 2 − 8w` 0 .
 
 
 
(2) (2)
0 2 + 5w` 1 − 3w`
Obtenemos que la matriz L = [L (α1 , β1 ) , L (α2 , β2 )] es
 
0 0 0 1 0 0
L= 0 0 0 0 0 0 .
0 0 1 0 0 0
La condición (2.228) se satisface para cualquier función vectorial {f (m)} de la forma
f (m) = (f1 (m), 0, f3 (m))T .

n o n o
Las proyecciones fb1 (m) , fb2 (m) definidas por (2.231) toman la forma
 
0
fb1 (m) = [L (α1 , β1 ) , 0] L† f (m) =  0  , (2.245)
f3 (m)
84 2.2 Caso general
 
f1 (m)
fb2 (m) = [0, L (α2 , β2 )] L† f (m) =  0  . (2.246)
0
Notemos que
1
 
1 − 0 

 2 
A + AH A + AH


 1
 1079 297
= − , σ = , , 1 ;

2  2 2 1   2 396 1079
 
 
0 1 1
 
−3 1 0
 
 
B + BH 
 1 −8 5 
B + BH

1211 1729

= , σ = − , −3, − ,

2  2 
 2 132 947
 
 5 
0 −3
2
y por tanto las condiciones de estabilidad (2.219), (2.220) se satisfacen. Tomando r sufi-
cientemente pequeño de forma que la condición (2.238) se verifica, y M verificando (2.168),
(2.184) y (2.238), por el teorema 2.7 la función vectorial
2
X
U (m, n) = Uδ (m, n) ,
δ=1
n o
donde {Uδ (m, n)} están definidas en (2.235), (2.236) y fbδ (m) por (2.245)–(2.246), es una
solución estable del problema mixto discreto (2.153)–(2.156).
Capı́tulo 3
Sistemas hiperbólicos
0.4pt0.4pt 0pt0.4pt
3.1. Introducción
En este capı́tulo utilizaremos esquemas matriciales en diferencias finitas para construir
soluciones numéricas discretas de problemas mixtos de tipo hiperbólico modelizados por
Auxx (x, t) − utt (x, t) = 0, 0 < x < 1, t > 0, (3.1)

u(0, t) = 0, t > 0, (3.2)
Bu(1, t) + Cux (1, t) = 0, t > 0, (3.3)
u(x, 0) = F (x), 0 ≤ x ≤ 1, (3.4)
ut (x, 0) = V (x), 0 ≤ x ≤ 1. (3.5)
donde A, B, C son s × s matrices complejas, elementos de Cs×s , y la incógnita u y F , V ,

son vectores s-dimensionales, elementos de Cs . Asumiremos que
C invertible , (3.6)

Es importante notar que incluso en el caso donde A es una matriz diagonalizable el pro-
blema permanece acoplado si las matrices B y C no son simultáneamente diagonalizables
con A.
En el apartado 3.2 estudiaremos la discretización del problema (3.1)–(3.5) utilizando apro-
ximaciones en diferencias centradas para las derivadas de segundo orden uxx , utt , aproxima-
ciones en diferencias progresivas para ut y en diferencias regresivas para ux . En el apartado
3.3 estudiaremos la existencia y construcción de soluciones no triviales para el problema
de contorno sin considerar las condiciones iniciales, utilizando un método de separación de
variables discreto. En el apartado 3.4 construiremos soluciones no triviales para el problema
mixto discreto haciendo uso del método de separación de variables discreto del apartado
3.3. En el apartado 3.5 comprobaremos que el esquema en diferencias finitas utilizado en el
apartado 3.2 es consistente respecto a la ecuación (3.1) y analizaremos la estabilidad de las
85
86 Capı́tulo 3. Sistemas hiperbólicos
soluciones construidas para el problema mixto discreto. En el apartado 3.6 utilizaremos un

método de proyecciones que nos permitirá extender los resultados del apartado 3.4 a clases
más amplias de funciones de valores iniciales F (x) y V (x). Finalmante, en el apartado 3.7
desarrollaremos algunos ejemplos para ilustrar los resultados teóricos obtenidos.
3.2. La discretización del problema

Dividimos el dominio [0, 1] × [0, ∞[ en rectángulos iguales de lados ∆x = h y ∆t = k,
introducimos las coordenadas de un punto tı́pico de la malla (mh, nk) y representamos
U (m, n) = u(mh, nk). Utilizando aproximaciones en diferencias centradas para uxx y utt ,
diferencias progresivas para ut y regresivas para ux (véase el apartado 2.1.2):
U (m, n + 1) − U (m, n)
ut (mh, nk) ≈ ,
k
U (m, n + 1) − 2U (m, n) + U (m, n − 1)
utt (mh, nk) ≈ ,
k2
U (m, n) − U (m − 1, n)
ux (mh, nk) ≈ ,
h
U (m + 1, n) − 2U (m, n) + U (m − 1, n)
uxx (mh, nk) ≈ ,
h2
obtenemos la discretización del problema (3.1)–(3.5):
r2 A [U (m + 1, n) + U (m − 1, n)] + 2(I − r2 A)U (m, n)

− [U (m, n + 1) + U (m, n − 1)] = 0 ,
1 ≤ m ≤ M − 1, n > 0 (3.7)
U (0, n) = 0, n>0 (3.8)

BU (M, n) + M C [U (M, n) − U (M − 1, n)] = 0, n>0 (3.9)
U (m, 0) = F (mh) = f (m) , 0 ≤ m ≤ M (3.10)
U (m, 1) − U (m, 0)
= V (mh) = v(m) , 0 ≤ m ≤ M (3.11)
k
donde
k 1
r= , h= . (3.12)
h M
3.3. El problema de contorno discreto

En este apartado buscamos soluciones no triviales del problema de contorno (3.7)–(3.9) de
la forma:
Capı́tulo 3. Sistemas hiperbólicos 87
U (m, n) = G(n)H(m) , G(n) ∈ Cs×s , H(m) ∈ Cs . (3.13)

Exigiendo la condición (3.13), el esquema en diferencias finitas dado en (3.7) se escribe de
la forma
r2 AG(n) [H(m + 1) + H(m − 1)] + 2 I − r2 A G(n)H(m)

− [G(n + 1) + G(n − 1)] H(m) = 0 . (3.14)
Tomamos ρ un número real y escribimos la ecuación (3.14) de la siguiente forma
h ρ i
r2 AG(n) H(m + 1) − 2 + 2 H(m) + H(m − 1)
r
− [G(n + 1) − (2I + ρA)G(n) + G(n − 1)] H(m) = 0 . (3.15)
Notemos que la ecuación (3.15) se satisface si las sucesiones {G(n)} y {H(m)} cumplen:
ρ
H(m + 1) − 2 + 2 H(m) + H(m − 1) = 0 , 1 ≤ m ≤ M − 1 , (3.16)
r
G(n + 1) − (2I + ρA) G(n) + G(n − 1) = 0 , n > 0 . (3.17)
Tomemos ρ ∈ R tal que
−4r2 < ρ < 0 . (3.18)

2
2r 2 +ρ
Entonces 2r2
< 1 y existe θ ∈]0, π[ tal que
2r2 + ρ

2 2 θ 2
cos θ = , ρ = 2r (cos θ − 1) = −4r sen . (3.19)
2r2 2
Sean
s
2r2 + ρ 2r + ρ 2
2
z0 = +i 1− = eiθ ,
2r2 2r2
s
2 2r + ρ 2
2
2r + ρ
z1 = − i 1 − = e−iθ ,
2r2 2r2
con i2 = −1, las soluciones de la ecuación escalar

ρ
z2 − 2 + 2 z + 1 = 0 , (3.20)
r
y observemos que
z0n = cos(nθ) + i sen(nθ) ; z1n = cos(nθ) − i sen(nθ) . (3.21)
La solución general de la ecuación (3.16) verificando que H(0) = 0 es de la forma
H(m) = sen(mθ)E , E ∈ Cs , 1 ≤ m ≤ M − 1. (3.22)

Por tanto se verifica la condición de contorno (3.8) y la condicón (3.9) se satisface si
L(θ)G(n)E = 0 , n > 0, (3.23)

donde {G(n)} es una solución de (3.17) y L(θ) es la matriz definida por
L(θ) = B sen(M θ) + M C [sen(M θ) − sen ((M − 1)θ)] . (3.24)

Para resolver (3.17) consideramos la ecuación matricial algebraica
W 2 − (2I + ρA) W + I = 0 , W ∈ Cs×s . (3.25)

Suponiendo que
A tiene todos los valores propios positivos , (3.26)
tomando
1 n 1 o
r<p = mı́n a− 2 ; a ∈ σ(A) , (3.27)
ρ(A)
θ
y utilizando que ρ = −4r2 sen2

2 , se sigue que las matrices
s 2 
ρA ρA 
−I 

W0 = I + + I+ 
2 2 


s , (3.28)
2 

ρA ρA 

W1 = I + − I+ −I 

2 2

son soluciones de (3.25) tales que

s 2
ρA
W0 − W1 = 2 I+ −I es invertible, (3.29)
2
i.e., W0 y W1 forman un sistema fundamental de soluciones de (3.25).
En efecto, por el teorema de la aplicación espectral, véase el teorema 1.6, se tiene
( r )
ρa 2
σ (W0 − W1 ) = 2 1+ − 1 ; a ∈ σ(A) (3.30)
2
y
ρa 2 ρ2 a2 ρa
1+ −1= + ρa = ρa 1 + 6 0,
= (3.31)
2 4 4
ya que se verifican las condiciones (3.26) y (3.27).
Por (3.29), véasen [18], [22], la solución general de (3.17) viene dada por
G(n) = W0n P + W1n Q P, Q ∈ Cs×s . (3.32)

Por las propiedades del cálculo funcional matricial las matrices W0 y W1 son polinomios
en la matriz A de grado p − 1, véase el teorema 1.10, donde
p es el grado del polinomio minimal de A . (3.33)

Por tanto, la condición (3.23) es equivalente a la condición
L(θ)Aj (P, Q) = 0 , 0 ≤ j < p , P, Q ∈ Cs . (3.34)

Por tanto, el problema de valores frontera (3.7)–(3.9) admite soluciones no triviales de la
forma (3.13) si existen vectores P , Q no nulos simultáneamente, verificando (3.34). Una
condición suficiente para tener soluciones no nulas P y Q es que la matriz
L(θ) sea singular . (3.35)

`π
Obsérvese que si θ` = M para 1 ≤ ` ≤ M − 1, entonces sen(M θ) = 0 y se obtiene que
`π
en este caso la matriz L(θ` ) = M C(−1)` sen M

es invertible ya que se verifica (3.6) y
`π
sen( M ) 6= 0. Por tanto los valores de θ para los cuales se verifica (3.35) deben satisfacer
que el sen(M θ) 6= 0 y podemos escribir (3.35) de la forma
M [sen(M θ) − sen((M − 1)θ)]

C −1 B + I es singular . (3.36)
sen(M θ)
Asumamos que
Existe µ ∈ σ(−C −1 B) ∩ R . (3.37)
La condición (3.36) se cumple si existen θ para la ecuación escalar
M [sen(M θ) − sen ((M − 1)θ)]

= µ. (3.38)
sen(M θ)
La ecuación (3.38) es equivalente a

µ
cos θ − 1 − M
cot(M θ) = . (3.39)
sen θ
Por el apartado 2.1.3 del capı́tulo 2.1 sabemos que
i h
(`−1)π `π
• Si µ < 0, existe una solución θ` de la ecuación (3.39) en J` = M ,M para
1 ≤ ` ≤ M − 1, de forma que la matriz L(θ` ) es singular.
• Si µ = 0, entonces B es singular y la ecuación (3.38) es equivalente a la ecuación
sen(M θ) = sen ((M − 1)θ) ,

(2`−1)π
cuyas soluciones son θ` = 2M −1 , para 1 ≤ ` ≤ M − 1.
i h
• Si µ > 0, existe una solución θ` de la ecuación (3.39) en J` = (`−1)π
M , `π
M para
2 ≤ ` ≤ M − 1, de forma que la matriz L(θ` ) es singular. Por tanto en este caso sólo
podemos garantizar analı́ticamente la existencia de M − 2 soluciones de la ecuación
(3.39).
Definiendo la matriz
C −1 L(θ)
G(µ) = = C −1 B + µI , (3.40)
sen(M θ)
obtenemos entonces que la condición (3.34) es equivalente a
G(µ)Aj (P, Q) = 0 0 ≤ j < p. (3.41)

Si definimos la matriz por bloques G(µ)
e como
 
G(µ)

 G(µ)A 

G
e (µ) = 
 G(µ)A2 ,

(3.42)
 .. 
 . 
G(µ)Ap−1
entonces (3.41) puede escribirse de la forma
G(µ)(P,
e Q) = 0 . (3.43)
Por el teorema de Mitra,véase el teorema 1.1, el sistema algebraico (3.43) admite soluciones
no nulas (P, Q) si rango G(µ)
e < s, y en este caso el conjunto solución de (3.43) se expresa
por
e † G(µ)

(P, Q) = I − G(µ) e (P0 , Q0 ) , P0 , Q0 ∈ Cs .
Resumiendo, bajo las hipótesis (3.6) y (3.37), si θ` son las soluciones de (3.39) y G(µ)
e la
matriz definida por (3.42), el conjunto de soluciones no triviales del problema de contorno
(3.7)–(3.9) viene dado por
U` (m, n) = (W0n P` + W1n Q` ) sen(mθ` )




1 ≤ m ≤ M − 1, n > 0






r 

2 

2 2 θ` 2 2 θ`
W0 = I − 2Ar sen 2 + I − 2Ar sen 2 −I 



, (3.44)
r 

2 
θ` θ`
W1 = I − 2Ar2 sen2

2 − I − 2Ar2 sen2 2 − I 








e † G(µ)

(P` , Q` ) = I − G(µ) e (P, Q) , s
P, Q ∈ C 

donde 1 ≤ ` ≤ M − 1 si µ ≤ 0 y 2 ≤ ` ≤ M − 1 si µ > 0.
Supongamos ahora que el parámetro real introducido en la ecuación (3.15) toma el valor
ρ = 0. En este caso la ecuación caracterı́stica (3.20) toma la forma:
z 2 − 2z + 1 = (z − 1)2 = 0,
y el conjunto de soluciones de la ecuación vectorial (3.16) viene dado por:

H(m) = 1m e
c + m 1m de = e
c + m d,
e c, de ∈ Cs , 1 ≤ m ≤ M − 1.
e (3.45)
Teniendo en cuenta que es necesario que H(0) = 0, para no tener soluciones triviales del
problema de frontera discreto (3.7)–(3.9), se obtiene
H(m) = m d,
e de ∈ Cs , 1 ≤ m ≤ M. (3.46)
Por otro lado, puesto que ρ = 0, la ecuación matricial (3.17) toma la forma:
G(n + 1) − 2G(n) + G(n − 1) = 0 , n > 0 , G(n) ∈ Cs×s . (3.47)

Por tanto una solución para la ecuación (3.47) viene dada por
G(n) = I n D
e + nI n Fe = D
e + nFe, e Fe ∈ Cs×s n > 0.
D, (3.48)
Si sustituimos (3.46) y (3.48) en la condición de contorno (3.9) llegamos a la siguiente
expresión:
BG(n)H(M ) + M C [G(n)H(M ) − G(n)H(M − 1)] = 0 ,
o equivalentemente

(B + C) De + nFe de = 0 , e Fe ∈ Cs×s , de ∈ Cs .
B, C, D, (3.49)
Utilizando que C es invertible podemos expresar (3.49) de la forma

(C −1 B + I) (c + nd) = 0 , c, d ∈ Cs , n > 0 . (3.50)
En vista de (3.50), el problema de contorno discreto (3.7)–(3.9) admitirá soluciones no
triviales de la forma (3.13) si existen vectores c, d, no nulos simultáneamente verificando
(C −1 B + I) (c, d) = 0 c, d ∈ Cs . (3.51)
Una condición suficiente para que se verifique (3.51) es que la matriz G(1) = C −1 B + I sea
singular, o equivalentemente
Exista µ = 1 ∈ σ −C −1 B .

(3.52)
Bajo la hipótesis (3.52), la ecuación
G(1)(c, d) = 0, c, d ∈ Cs , (3.53)
tiene soluciones c y d no nulas si y sólo si rango (G(1)) < s, y bajo esta condición por el
teorema de Mitra, véase el teorema 1.1, el conjunto solución de (3.53) viene dado por

(c, d) = I − G(1)† G(1) (c0 , d0 ) , c0 , d0 ∈ Cs .
Por tanto, para µ = 1 obtenemos la solución no trivial del problema de contorno (3.7)–(3.9):

U (m, n) = m (c + n d) , 1 ≤ m ≤ M − 1, n > 0 
. (3.54)
(c, d) = I − G(1)† G(1) (c0 , d0 ) , c0 , d0 ∈ Cs ∼ {0}


TEOREMA 3.1 Consideremos el problema de frontera (3.7)–(3.9) y asumamos la con-

dición (3.6). Sea µ un número real verificando (3.37). Sean G(µ) y G e (µ) las matrices en
s×s k
C definidas por (3.40) y (3.42) respectivamente. Sea M > 0; r = h verificando la condi-
ción (3.27) y p el grado del polinomio minimal de la matriz A.
(i) Si µ ∈ R, existen valores θ ∈ ]0, π[ tales que la matriz L(θ) definida por (3.24) es
singular. Bajo la hipótesis (3.35) el problema
de contorno (3.7)–(3.9) tiene soluciones
no triviales de la forma (3.13) si rango G(µ)e < s.

2`−1
(ii) Si C es singular, se toma el valor propio µ = 0 entonces θ` = 2M −1 π para 1 ≤ ` ≤

−1
M − 1 son soluciones de (3.39). Si rango G(0)
e < s entonces {U` (m, n)}M
`=1 dadas
por (3.44) define M − 1 soluciones no triviales del problema de contorno (3.7)–(3.9).
(iii) Si µ < 0, entonces

i parah 1 ≤ ` ≤ M − 1 existe una solución θ` de la ecuación
(`−1)π `π
(3.39) en J` = M , M para la cual la matriz L(θ` ) es singular. Bajo la hipótesis

−1
rango G(µ)
e < s, la sucesión {U` (m, n)}M
`=1 dada por (3.44) define M −1 soluciones
no triviales del problema de contorno (3.7)–(3.9).
(iv) Si µ = 1, entonces
i parah 2 ≤ ` ≤ M − 1 existe una solución θ` de la ecuación
(`−1)π `π
(3.39) en J` = M , M para la cual la matriz L(θ` ) es singular. Bajo las hipótesis

rango G(1)
e < s y rango (G(1)) < s, existen dos conjuntos de soluciones no
triviales del problema de contorno (3.7)–(3.9) dadas por (3.44) y (3.54).
(v) Si µ ∈ ]0, +∞[∼ {1}, entonces

i h 2 ≤ ` ≤ M − 1 existe una solución θ` de la
para
(`−1)π `π
ecuación (3.39) en J` = M , M para la cual la matriz L(θ` ) es singular. Bajo la

−1
hipótesis rango G(µ)
e < s, la sucesión {U` (m, n)}M
`=2 dada por (3.44) define M − 2
soluciones no triviales del problema de contorno (3.7)–(3.9).
3.4. El problema mixto discreto

En este apartado construiremos una solución exacta del problema mixto (3.7)–(3.11). Asu-
miremos las hipótesis y la notación del apartado 3.3.
Por el apartado 2.1.5 del capı́tulo 2.1 sabemos que las condiciones de contorno (3.8)–(3.9)
y las hipótesis (3.6) y (3.37) nos conducen al siguiente problema de Sturm-Liouville escalar
discreto
ρ

−h(m + 1) + 2h(m) − h(m − 1) = − 2 h(m) , 

r . (3.55)
M
h(0) = 0 ; h(M ) = h(M − 1) , 1 ≤ m ≤ M − 1 
M −µ

exactamente M −1 valores
Además en el apartado 2.1.5, vimos que el problema (3.55) tiene
−ρ` M −1
propios reales dados por r2 `=1 , donde ρ` = −4r sen θ2` y θ` verifica (3.39). Para
2 2
cada valor propio −ρ

r2
`
existe una función propia {h` (m)} = {sen(mθ` )}. Estas funciones
propias son ortogonales con respecto a la función peso 1, véase el apartado 1.2.
Notemos que se ha considerado el problema de Sturm-Liouville discreto escalar (3.55) a

fin de identificar los vectores P` , Q` y c, d que aparecen en las soluciones (3.44) y (3.54)
respectivamente.
CASO 1. µ=0 (B singular)

2`−1
Supongamos que rango G(0) < s , θ` = 2M −1 π , 1 ≤ ` ≤ M − 1. Por superposición
e
de las soluciones del problema de contorno (3.7)–(3.9) obtenidas en el teorema 3.1-(ii), se
sigue que la sucesión vectorial
M
X −1
U (m, n) = (W0n P` + W1n Q` ) sen (mθ` ) , 1 ≤ m ≤ M − 1, n > 0, (3.56)
`=1
es una solución del problema de contorno (3.7)–(3.9), donde los vectores P` y Q` están en
Nuc G(µ)
e y deben ser elegidos de forma que se satisfagan las condiciones (3.10) y (3.11).
Imponiendo las condiciones iniciales (3.10) y (3.11) a la sucesión vectorial U (m, n) definida
por (3.56) se sigue que
M
X −1
f (m) = (P` + Q` ) sen (mθ` ) , (3.57)
`=1
M
X −1
kv(m) + f (m) = {W0 P` + W1 Q` } sen (mθ` ) . (3.58)
`=1
Por la teorı́a de las series de Fourier discretas, véase el apartado 1.2, trabajando componente
a componente en (3.57) y (3.58) se sigue que
MP−1
sen(νθ` )f (ν)
ν=1
P` + Q` = M −1
, (3.59)
sen2 (νθ
P
`)
ν=1
MP−1
{kv(ν) + f (ν)} sen(νθ` )
ν=1
W0 P` + W1 Q` = M −1
. (3.60)
sen2 (νθ
P
`)
ν=1
Premultiplicando (3.59) por W0 y restando (3.60) obtenemos

MP−1
{(W0 − I) f (ν) − kv(ν)} sen(νθ` )
ν=1
(W0 − W1 ) Q` = M −1
. (3.61)
sen2 (νθ
P
`)
ν=1
Premultiplicando (3.59) por (−W1 ) y sumando (3.60) se sigue que
MP−1
{kv(ν) − (W1 − I) f (ν)} sen(νθ` )
ν=1
(W0 − W1 ) P` = M −1
. (3.62)
sen2 (νθ
P
`)
ν=1
Puesto que W0 − W1 es invertible, de (3.61)–(3.62) obtenemos
M −1
(W0 − W1 )−1
P
{kv(ν) − (W1 − I) f (ν)} sen(νθ` )
ν=1
P` = M −1
, (3.63)
sen2 (νθ` )
P
ν=1
M −1
(W0 − W1 )−1
P
{(W0 − I) f (ν) − kv(ν)} sen(νθ` )
ν=1
Q` = M −1
. (3.64)
sen2 (νθ
P
`)
ν=1
Ya que por (3.28), (3.29) las matrices W0 , W1 y (W0 − W1 )−1 son polinomios en la matriz
A de grado p − 1 (véase el teorema 1.10), por (3.42), (3.63), (3.64) los vectores P` y Q`
satisfacen (3.43) si
{f (m), v(m) , 1 ≤ m ≤ M − 1} ⊂ Nuc G(0) , (3.65)
Nuc G(0) es un subespacio invariante por la matriz A , (3.66)
véase la definición 1.2. La condición (3.66) se puede expresar de la forma

G(0)A I − G(0)† G(0) = 0 , (3.67)

Por tanto, bajo las hipótesis (3.65) y (3.67), la expresión

M −1
(W0n P` W1n Q` ) sen (mθ` )
P
U (m, n) = + , 



`=1 





M −1


−1 P 
(W0 − W1 ) {kv(ν) − (W1 − I) f (ν)} sen(νθ` ) 


ν=1 
P` = M −1
, 



2
P
sen (νθ` )


ν=1 , (3.68)



M −1


−1 P 
(W0 − W1 ) {(W0 − I) f (ν) − kv(ν)} sen(νθ` ) 

ν=1


Q` = M −1
, 



2
P
sen (νθ` )




ν=1  

2`−1
θ` = 2M −1 π 

1 ≤ m ≤ M − 1, n > 0, 1 ≤ ` ≤ M − 1
donde W0 y W1 están definidas en (3.44), es una solución del problema mixto discreto
(3.7)–(3.11).
CASO 2. µ<0

−1
Sea µ ∈ σ(−C −1 B) y supongamos que rango G(µ)
e < s donde {θ` }M
`=1 vienen dadas por
el teorema 3.1-(iii). Bajo las hipótesis
{f (m), v(m) , 1 ≤ m ≤ M − 1} ⊂ Nuc G(µ) , (3.69)
Nuc G(µ) es subespacio invariante por la matriz A , (3.70)
o equivalentemente

G(µ)A I − G(µ)† G(µ) = 0 , (3.71)
se obtiene que la expresión


M −1
P
U (m, n) = + , 



`=1 





M −1


−1 P 
(W0 − W1 ) {kv(ν) − (W1 − I) f (ν)} sen(νθ` ) 


ν=1 
P` = M −1
, 


2
P
sen (νθ` ) , (3.72)
ν=1 





M −1


−1 P 
(W0 − W1 ) {(W0 − I) f (ν) − kv(ν)} sen(νθ` ) 

ν=1


Q` = M −1
, 



2
P
sen (νθ` )



ν=1
1 ≤ m ≤ M − 1, n > 0, 1 ≤ ` ≤ M − 1
donde W0 y W1 están definidas en (3.44), es una solución del problema mixto discreto
(3.7)–(3.11).
CASO 3. µ=1

Supongamos que rango G(1)
e < s, y rango (G(1)) < s. Por superposición de las soluciones
del problema de contorno (3.7)–(3.9) obtenidas en el teorema 3.1-(iv), se sigue que la
sucesión vectorial
M −1

P
U (m, n)= + + m (c + n d) ,
`=2 (3.73)
(P` , Q` ) ∈ Nuc G(1)
e , (c, d) ∈ Nuc G(1)
1 ≤ m ≤ M − 1, n > 0, 2 ≤ ` ≤ M − 1
es una solución del problema de contorno (3.7)–(3.9). Imponiendo las condiciones iniciales
(3.10) y (3.11) a la solución (3.73) obtenemos que los vectores (P` , Q` ) y (c, d) deben
verificar
M
X −1
f (m) = (P` + Q` ) sen (mθ` ) + mc , (3.74)
`=2
M
X −1
kv(m) + f (m) = (W0 P` + W1 Q` ) sen (mθ` ) + m(c + d) . (3.75)
`=2
En el apartado 2.1.5, se probó que las funciones propias del problema escalar discreto de
Sturm-Liouville (3.55):
−1
h` (m) = {sen (mθ` )}M
`=2
M −1 n oM −1
asociadas a los valores propios − rρ2` `=2 = 4 sen2 θ2`

,y
`=2
h1 (m) = m
asociada al valor propio 0, son ortogonales respecto a la función peso 1. Por tanto en vista
de (3.74) y razonando como en el apartado 2.1.5 obtenemos
MP−1
sen(νθ` )f (ν)
ν=1
P` + Q` = M −1
, (3.76)
sen2 (νθ` )
P
ν=1
y
M −1
1 X
c= νf (ν) . (3.77)
M3 M2 M
3 − 2 + 6 ν=1
De la misma forma, de (3.75) se obtiene

MP−1
{kv(ν) + f (ν)} sen(νθ` )
ν=1
W0 P` + W1 Q` = M −1
, (3.78)
sen2 (νθ
P
`)
ν=1
y
M −1
1 X
c+d= ν {kv(ν) + f (ν)} . (3.79)
M3 M2 M
3 − 2 + 6 ν=1
Razonando análogamente al Caso 1 de este apartado, de (3.76) y (3.78) obtenemos los

−1
vectores {P` , Q` }M
`=2 definidos en (3.63) y (3.64) respectivamente. Sustituyendo (3.77) en
(3.79) obtenemos
M −1
k X
d= 3 ν v(ν) . (3.80)
M M2 M
3 − 2 + 6 ν=1
Además como se tiene que verificar que (P` , Q` ) ∈ Nuc G(1) e y (c, d) ∈ Nuc G(1), en vista
de (3.63), (3.64), (3.77) y (3.80) será suficiente exigir las condiciones (3.69) y (3.71).
CASO 4. µ ∈ ]0, +∞[, µ 6= 1
En este caso razonado

análogamente al Caso 4 del apartado 2.1.5 obtenemos, bajo las con-
diciones rango G(1) < s, (3.69) y (3.71), las M − 1 soluciones dadas por (3.72).
e
TEOREMA 3.2 Bajo las hipótesis (3.69) y (3.71) junto con el teorema 3.1, se obtiene:
(i) Si µ ∈ R ∼ {1}, la sucesión vectorial {U (m, n)} definida por (3.72) es una solución
exacta del problema mixto (3.7)–(3.11), donde para el caso µ = 0 se tiene θ` =

2`−1
2M −1 π, 1 ≤ ` ≤ M − 1.
(ii) Si µ = 1, la sucesión vectorial {U (m, n)} definida por (3.73), (3.63), (3.64), (3.77)
y (3.80) es una solución exacta del problema mixto (3.7)–(3.11).
3.5. Consistencia y estabilidad de las soluciones

Veamos que el esquema en diferencias (3.7) es consistente respecto a la ecuación (3.1) en
el sentido de la definición 1.13. Denotemos:
Λ[u] = utt (x, t) − Auxx (x, t) , (3.81)
u (mh, (n + 1)k) − u(mh, nk) + u(mh, (n − 1)k)

Λk,h [u] =
k2 (3.82)
[u ((m + 1)h, nk) − 2u(mh, nk) + u ((m − 1)h, nk)]
−A .
h2
Sea φ(x, t) una función suficientemente derivable y sea Φ(m, n) = φ(mh, nk). Considerando
el desarrollo de Taylor de Φ alrededor del punto (m, n) para: (m, n + 1), (m, n − 1) con
respecto a la variable t y para los puntos: (m + 1, n) y (m − 1, n) con respecto a la variable
x, se sigue que:
2 3 
Φ(m, n+1) = Φ(m, n)+kΦt(m, n)+ k2! Φtt(m, n)+ k3! Φttt(m, n)+ O(k 4 )

2 3

Φ(m, n−1) = Φ(m, n)−kΦt(m, n)+ k2! Φtt(m, n)− k3! Φttt(m, n)+
O(k )   4
(3.83)


h2 h3
Φ(m+1, n) = Φ(m, n)+hΦx(m, n)+ 2! Φxx(m, n)+ 3! Φxxx(m, n)+ O(h4 )



2 3
h h 4
Φ(m−1, n) = Φ(m, n)−hΦx(m, n)+ 2! Φxx(m, n)− 3! Φxxx(m, n)+ O(h )
Sustituyendo (3.83) en Λk,h [φ], dada por (3.82), obtenemos:
Λk,h [φ] = Φtt (m, n) − AΦxx (m, n) + O(k 2 ) + O(h2 ) . (3.84)
Teniendo en cuenta la expresión de Λ[φ], dada por (3.81), y la expresión (3.84) se obtiene:
kΛ[φ] − Λk,h [φ]k ≤ O(k 2 ) + O(h2 ) −→ 0 , cuando k, h → 0 . (3.85)
Por tanto, el esquema en diferencias (3.7) es consistente con la ecuación (3.1).

Ahora estudiaremos la estabilidad de las soluciones discretas construidas en el teorema 3.2.

Por el teorema de la aplicación espectral, véase el teorema 1.6, el espectro de las matrices
W0 y W1 definidas en (3.44) viene dado por
 
s 2
 θ` θ` 
σ (W0 ) = 1 − 2ar2 sen2 + 1 − 2ar2 sen2 − 1 : a ∈ σ(A) ,
 2 2 
y  
s 2
 θ` θ` 
σ (W1 ) = 1 − 2ar2 sen2 − 1 − 2ar2 sen2 − 1 : a ∈ σ(A) .
 2 2 
No podemos garantizar la estabilidad en el sentido de la definción 1.14, es decir, no podemos

exigir que las matrices W0 y W1 sean convergentes, véase el lema 1.1, i.e. no podemos exigir
que:
ρ (W0 ) , ρ (W1 ) < 1 ,
o equivalentemente

s 2
2
1 − 2ar sen2 θ ` 2 2
θ`
± 1 − 2ar sen −1 < 1,

2 2

ya que

s 2
1 − 2ar2 sen2 θ` + θ`

2
1 − 2ar sen 2 − 1 ·

2 2
s
2
2 2 θ` θ`

· 1 − 2ar sen
− 2
1 − 2ar sen 2 − 1 = 1 .
2 2
Por tanto analizaremos la estabilidad de las soluciones construidas en el teorema 3.2 en el

sentido de la definición 1.15. Esto significa que dado un punto del mallado (X, T ) donde
m 1
X = M = mh0 , h0 = M fijo, y T = Jk finito, estamos interesados en el comportamiento
de las soluciones {U (m, n)} cuando k → 0, i.e., cuando n crece, dependiendo del valor de
k, hasta tomar el valor J necesario para alcanzar el tiempo T . Por (3.26) y por la forma
de las matrices W0 y W1 definidas en (3.44) se sigue

kW0 k ≤ 1 + O(r); kW0 − Ik = O(r); (W0 − W1 )−1 = O(r−1 ),


(3.86)

kW1 k ≤ 1 + O(r); kW1 − Ik = O(r); r → 0

1
Puesto que h0 = M está fijado y r = hk , la expresión (3.86) significa

kW0 k ≤ 1 + O(k); kW0 − Ik = O(k); (W0 − W1 )−1 = O(k −1 ),


(3.87)

kW1 k ≤ 1 + O(k); kW1 − Ik = O(k); k → 0

Por (3.63), (3.64) y (3.87) se sigue que
kP` k = O(1), kQ` k = O(1), k → 0 . (3.88)

Por (3.88) se obtiene que {U (m, n)} permanece acotada cuando n crece, si los números
kW0n k , kW1n k permanecen acotados cuando n → ∞

.
siendo 0 ≤ n ≤ J, T = Jk finito
En efecto, ya que se tiene kW0 k ≤ 1 + O(k) podemos tomar kW0 k ≤ 1 + kS, para alguna
constante positiva S, y entonces para 0 ≤ n ≤ J obtenemos
kW0n k ≤ kW0 kn ≤ (1 + O(k))n ≤ (1 + O(k))J

≤ eJ O(k) ≤ eJkS = eT S = O(1). (3.89)
Lo mismo ocurre para kW1n k. Por tanto, por (3.88)–(3.89), las soluciones definidas en (3.72)
para µ ∈ R ∼ {1} son estables, i.e.,
1

kU (m, n)k = O(1) , 1 ≤ m ≤ M − 1, h0 = M fijo
. (3.90)
n → ∞, siendo 0 ≤ n ≤ J, T = Jk finito
Por otro lado para el caso µ = 1 también tenemos que la solución {U (m, n)} definida
por (3.73), (3.63), (3.64), (3.77) y (3.80) es estable. En efecto, en vista de la definición de
estabilidad, en el sentido de la definición 1.15, tenemos
km(c + nd)k ≤ k1 m, para alguna constante positiva k1 , 1 ≤ m ≤ M, (3.91)

puesto que se considera T = Jk finito con 0 ≤ n ≤ J, y los vectores (c, d) verifican
kck = O(1) , kdk = O(k) , k → 0.

1
Además como se fija el paso espacial h0 = M y 1 ≤ m ≤ M , de (3.91) podemos concluir
que
km(c + nd)k = O(1) , k → 0. (3.92)

TEOREMA 3.3 Bajo las hipótesis del teorema 3.2, tomando los vectores f (m) y v(m)
acotados, y la matriz A y el parámetro r > 0 verificando (3.26) y (3.27) respectivamente
se tiene:
(i) Si µ ∈ R ∼ {1}, la sucesión vectorial {U (m, n)} definida por (3.72) es una solución
estable del problema mixto (3.7)–(3.11), en el sentido de (3.90).
(ii) Si µ = 1, la sucesión vectorial {U (m, n)} definida por (3.73), (3.63), (3.64), (3.77)
y (3.80) es una solución estable del problema mixto (3.7)–(3.11), en el sentido de
(3.90)–(3.92).
3.6. El método de las proyecciones

En este apartado abordaremos la construcción de soluciones del problema mixto discreto
(3.7)–(3.9) para funciones f (m) y v(m) satisfaciendo condiciones más generales que las
exigidas en el apartado 3.4. Supongamos que existen
{µ1 , . . . , µq } ⊂ σ −C −1 B ∩ R ,

(3.93)
donde µj 6= µk para 1 ≤ j, k ≤ q, j 6= k. Introducimos la matriz G (µj ) en Cs×s
G (µj ) = C −1 B + µj I , 1 ≤ j ≤ q. (3.94)
Como los polinomios x − µj son mútuamente coprimos (primos 2 a 2), véase la defini-
ción 1.11, entonces por el teorema de descomposición, véase el teorema 1.9, si L(x) es el
polinomio
L(x) = (x − µj ) (x − µ1 ) (x − µ2 ) · · · (x − µj−1 ) (x − µj+1 ) · · · (x − µq ) , (3.95)

entonces
S = Nuc L −C −1 B = Nuc G (µ1 ) ⊕ Nuc G (µ2 ) ⊕ · · · ⊕ Nuc G (µq ) .

(3.96)
Consideramos la sucesión de polinomios coprimos de grado q − 1 definidos por

q
Y
Qj (x) = (x − µk ) , 1≤j≤q. (3.97)
k=1,k6=j
Entonces por el teorema de Bezout, véase el teorema 1.8, tomando los escalares αj de la
forma
 −1
q
 Y 
αj = (µj − µk ) , 1≤j≤q,
 
k=1,k6=j
obtenemos
q
X
1 = Q(x) = αj Qj (x) . (3.98)
j=1
Considerando el polinomio
T (µj ) = G (µ1 ) · · · G (µj−1 ) G (µj+1 ) · · · G (µq ) , (3.99)
y aplicando el cálculo funcional matricial sobre la matriz −C −1 B, de (3.94), (3.95), (3.98)

y (3.99) se sigue que
q
X
I = (−1)q−1 αj T (µj ) . (3.100)
j=1
Por tanto consideraremos las proyecciones de f (m) y v(m) sobre el subespacio Nuc G (µj )
definidas por

fbj (m) = (−1)q−1 αj T (µj )f (m) 
, 1 ≤ j ≤ q, 0 ≤ m ≤ M (3.101)
vbj (m) = (−1)q−1 αj T (µj )v(m)

donde por (3.100) y (3.101) se satisface

q
X q
X
f (m) = fbj (m) , v(m) = vbj (m) , 0 ≤ m ≤ M . (3.102)
j=1 j=1
Supongamos que f (m), v(m) ∈ S, 0 ≤ m ≤ M , i.e.,
L −C −1 B f (m) = L −C −1 B v(m) = 0 ,

0≤m≤M, (3.103)
o equivalentemente por (3.95), (3.101) y (3.103) que se cumple

n o
fbj (m) , vbj (m), 1 ≤ j ≤ q , 0 ≤ m ≤ M ⊂ Nuc G (µj ) . (3.104)
Por tanto consideraremos q problemas diferentes obtenidos al reemplazar las condiciones

iniciales (3.10)–(3.11) por:

U (m, 0) = fbj (m) 


, 1 ≤ j ≤ q, 0 ≤ m ≤ M . (3.105)
U (m, 1) − U (m, 0) 
= vbj (m) 

k
Sea (Pj ), 1 ≤ j ≤ q, el problema definido por (3.7)–(3.9) y (3.105). Observar que por resul-
tados previos, el cambio de las condiciones iniciales en cada problema (Pj ) sólo modifica a
los vectores (P` , Q` ) y (c, d) que aparecen en (3.63), (3.64), (3.77) ny (3.80),
o porque ahora
los coeficientes de Fourier discretos involucran a las proyecciones fbj (m) y {b vj (m)} en
lugar de a los vectores de las condiciones iniciales {f (m)} y {v(m)}. Notemos que para los
problemas (Pj ), los coeficientes de Fourier toman la forma
si µj ∈ R ∼ {1}, para 1 ≤ ` ≤ M − 1 se obtiene:


M −1 n o
−1 P (j)
(W0,j − W1,j ) vj (ν) − (W1,j − I) fbj (ν) sen(νθ` ) 
kb
(j) ν=1
P` =

M −1 
(j) 
sen2 νθ`
P

ν=1


,

M −1 n o
(j)

−1 P
(W0,j − W1,j ) (W0,j − I) fbj (ν) − kb
vj (ν) sen(νθ` ) 

(j) ν=1
Q` =

M −1 
(j)
sen2 νθ`
P 
ν=1
(3.106)
o si µj = 1, para 2 ≤ ` ≤ M − 1 se obtiene:
M −1 n o 
(j)
(W0,j − W1,j )−1
P
vj (ν) − (W1,j − I) fbj (ν) sen(νθ` ) 
kb
(j) ν=1
P` = M −1


(j)
sen2 νθ`
P 

ν=1 


M −1 n o 
−1 P (j) 
(W0,j − W1,j ) (W0,j − I) fbj (ν) − kb
vj (ν) sen(νθ` )  ,
(j) ν=1
Q` = M −1


(j)
sen2 νθ`
P 

ν=1 


M −1 M −1

c = M3 M2 M
(j) 1
ν fbj (ν) ; d = M 3 M 2 M
(j) k
P P 
ν vbj (ν)
3
− 2
+ 6 ν=1 3
− 2
+ 6 ν=1
(3.107)
donde el superı́ndice (j) y W0,j , W1,j denota que se está trabajando con el valor propio
µj . Denotando Uj (m, n) a la solución del problema mixto discreto (3.7)–(3.9) y (3.105)
asociado al valor propio µj , 1 ≤ j ≤ q, se obtiene que por construcción de las proyecciones
la solución del problema mixto discreto original (3.7)–(3.11) es de la forma:
q
X
U (m, n) = Uj (m, n) .
j=1

TEOREMA 3.4 Con la notación del teorema 3.1, sean µ1 , . . . , µq valores propios reales
y distintos de la matriz −C −1 B. Sean G(µj ) y T (µj ) matrices en Cs×s definidas por (3.94)
e j ) la matriz en Cps×s definida por
y (3.99). Sea G(µ
 
G(µj )
 G(µj )A 
Ge (µj ) =  , (3.108)

 ..
 . 
G(µj )Ap−1
siendo p el grado del polinomio minimal de A. Supongamos que las proyecciones fbj (m) y
vbj (m), definidas en (3.101), verifican (3.104); que la matriz A satisface la condición (3.26)
y además

G (µj ) A I − G (µj )† G (µj ) = 0 , 1 ≤ j ≤ q , (3.109)
y que el parámetro r satisface (3.27).
i h
(j) (`−1)π `π
(i) Si µj ∈ R ∼ {1}, se tiene que θ` ∈ J` = M ,M , 1 ≤ ` ≤ M − 1, denotan las
M − 1 soluciones de la ecuación

(j) µj
cos θ` − 1− M
(j)
cot M θ` = ,1 ≤ ` ≤ M − 1, (3.110)
(j)
sen θ`

(j) 2` − 1
donde θ` = π si µ = 0. Para este caso asumiremos que rango Ge (µj ) <
2M − 1
s.
i h
(j) (`−1)π `π
(ii) Si µj = 1, se tiene que θ` ∈ J` = M ,M , 2 ≤ ` ≤ M − 1, denotan las M − 2
soluciones de la ecuación

(j) µ
cos θ` − 1 − Mj
(j)
cot M θ` = ,2 ≤ ` ≤ M − 1, (3.111)
(j)
sen θ`

Para este caso asumiremos que rango G(1)
e < s y rango (G(1)) < s.
Entonces bajo todas las hipótesis consideradas anteriormente se obtiene que
q
X
U (m, n) = Uj (m, n) , 1 ≤ m ≤ M − 1, n > 0, (3.112)
j=1
es una solución convergente (por ser consistente y estable) a la solución del problema mixto
continuo (3.1)–(3.5), donde
si µj ∈ R ∼ {1}:
M −1 n o
(j) (j) (j)
X
n n
Uj (m, n) = W0,j P` + W1,j Q` sen m θ` ,
`=1
1 ≤ m ≤ M − 1, n > 0
siendo
s 2 
(j) (j)
θ` θ`
W0,j = I − 2Ar2 sen2 2 + I − 2Ar2 sen2
−I  2


s

 (3.113)
(j) (j) 2 
2 2 θ` 2 2 θ`
= I − 2Ar sen − I − 2Ar sen −I

W1,j 2 2
(j) M −1
n o
(j)
y P` , Q` los vectores definidos en (3.106);
`=1
o si µj = 1:
M −1 n o
(j) (j) (j)
X
n n
Uj (m, n) = W0,j P` + W1,j Q` sen m θ` + m c(j) + nd(j) ,
`=1
1 ≤ m ≤ M − 1, n > 0
(j) M −1
n o n o
n ,Wn (j)
siendo W0,j 1,j las matrices definidas en (3.113) para 2 ≤ ` ≤ M −1, y P` , Q`
`=2
y c(j) , d(j) los vectores definidos en (3.107).

3.7. Ejemplos
     
2 −1 1 −1 −1 1 1 1 −1
A= 1 1 0 , B =  1 1 −2 , C =  1 −1 0 ,
1 −1 2 1 0 −1 1 0 −1
f (m) = (f1 (m), f2 (m), f3 (m))T , v(m) = (v1 (m), v2 (m), v3 (m))T ∈ C3 ,
1
h= M , 1 ≤ m ≤ M − 1.
La matriz C es invertible y −C −1 B viene dada por

 
1 0 0
−C −1 B =  2 1 −2  .
2 0 −1
En este caso tenemos σ −C −1 B = {−1, 1}. Denotamos µ1 = −1, µ2 = 1. Veamos prime-

ramente si se verifican las siguientes condiciones

rango G e (µj ) < 3 , 1 ≤ j ≤ 2 , y rango (G(1)) < 3 .
e (µj ), 1 ≤ j ≤ 2
Puesto que el grado del polinomio minimal de A es p = 3 las matrices G
toman la forma
 
G (µj )
e (µj ) =  G (µj ) A  , 1 ≤ j ≤ 2 ,
G
G (µj ) A2
siendo
   
−2 0 0 0 0 0
G (µ1 ) = C −1 B − I =  −2 −2 2 , G(1) = C −1 B + I =  −2 0 2  .
−2 0 0 −2 0 2
Únicamente podemos trabajar con elvalor propio µ2 = 1, ya que para el valor propio
µ1 = −1 se verifica que rango G e (µ1 ) = 3. Estamos en el caso (ii) del teorema 3.3. Para
µ2 = 1 se satisface que el Nuc G(1) es subespacio invariante por la matriz A, es decir:

G(1)A I − G(1)† G(1) = 0 ,
siendo
0 − 81 − 18
 
 
†
 
 0
G(1) =  0 0 
,
 
1 1
0 8 8
y además σ(A) = {1, 2}, i.e., todos los valores propios de A son positivos.
Por otro lado, para que los vectores {f (m)} y {v(m)} satisfagan
{f (m), v(m) , 1 ≤ m ≤ M − 1} ⊂ Nuc G(1) ,
debemos exigir que sean de la forma
f (m) = (f1 (m), f2 (m), f1 (m))T , v(m) = (v1 (m), v2 (m), v1 (m))T .
Entonces tomando 0 < r < √ 1 = 21 , por (ii) del teorema 3.3 la solución numérica estable
ρ(A)
del problema viene dada por
M
X −1
U (m, n)= (W0n P` + W1n Q` ) sen (mθ` ) + m (c + n d) ,
`=2
M −1
donde los vectores {P` , Q` }`=2 y (c, d) están definidos en (3.63)–(3.64), (3.77) y (3.80)
respectivamente.

1 
−1 21 2
    
3 0 3 −3 1 −1
     
   2 7 
  
A =  −24 1 −4 , B =  − 3 1 3 , C =  1 2
    3 
,
     
1
6 0 2 0 2 1 0 1 2
f (m) = (f1 (m), f2 (m), f3 (m))T , v(m) = (v1 (m), v2 (m), v3 (m))T ∈ C3 ,
1
h= M , 1 ≤ m ≤ M − 1.
La matriz C es invertible y −C −1 B viene dada por

 1 
6 0 0
 
−1 1 5 
 
−C B =  1 − 2 − 3  .

 
1 1
−2 0 3
En este caso tenemos σ −C −1 B = − 12 , 16 , 13 . Denotamos µ1 = − 12 , µ2 = 1 1

6, µ3 = 3.
Veamos primeramente si se verifican las siguientes condiciones:

rango G e (µj ) < 3 , 1 ≤ j ≤ 3 ,
e (µj ), 1 ≤ j ≤ 3
Puesto que el grado del polinomio minimal de A es p = 2 las matrices G
toman la forma

G (µj )
G (µj ) =
e , 1 ≤ j ≤ 3,
G (µj ) A
siendo
− 23 1
     
0 0 0 0 0 6 0 0
     
5 1 5 5 5
     
G (µ1 ) =  −1 0

3 , G (µ2 ) =  −1

3 3
, G (µ3 ) =  −1
  6 3
.

     
1 5 1
2 0 −6 2 0 − 16 1
2 0 0
Únicamente podemos trabajar con los valores

propios
µ1 = − 12 , y µ2 = 16 , ya que para el
valor propio µ3 = 13 se verifica que rango Ge (µ3 ) = 3. Estamos en el caso del teorema 3.4.
Para µ1 = − 21 , y µ2 = 61 se satisface que el Nuc G (µj ), 1 ≤ j ≤ 2, es subespacio invariante
por la matriz A, es decir:

G (µj ) A I − G (µj )† G (µj ) = 0 , 1 ≤ j ≤ 2 ,
siendo
− 32 18 27
   
0 0 0 77 11
   
† † 15 6
   
G (µ1 ) = 
 0 0 0  , G (µ2 ) = 

 0 77 11
,

   
9 12 6 54 15
− 10 25 − 25 0 77 11
y además σ(A) = {1, 4}, i.e., todos los valores propios de A son positivos.
n o2
Por otro lado, para que las proyecciones fbj (m) vj (m)}2 cumplan
y {b j=1
j=1
n o
fbj (m), vbj (m) , 1 ≤ m ≤ M − 1 ⊂ Nuc G (µj ) , 1 ≤ j ≤ 2,
debemos exigir que sean de la forma

)
fb1 (m) = (0, f1,2 (m), 0)T
T ,
fb2 (m) = 31 f2,3 (m), −4f2,3 (m), f2,3 (m)
y
)
vb1 (m) = (0, v1,2 (m), 0)T
T ,
vb2 (m) = 31 v2,3 (m), −4v2,3 (m), v2,3 (m)
donde se ha denotado a los vectores fbj (m) y vbj (m) de la siguiente forma
fbj (m) = (fj,1 (m), fj,2 (m), fj,3 (m))T y vbj (m) = (vj,1 (m), vj,2 (m), vj,3 (m))T .
Entonces tomando 0 < r < 21 , por el teorema 3.4 la solución numérica estable del problema
viene dada por
2
X
U (m, n) = Uj (m, n) , 1 ≤ m ≤ M − 1, n > 0,
j=1
siendo
M −1 n o
(j) (j) (j)
X
n n
Uj (m, n) = W0,j P` + W1,j Q` sen m θ` ,
`=1
1 ≤ m ≤ M − 1, n > 0
(j) M −1
n o
(j)
W0,j y W1,j las matrices definidas en (3.113), P` , Q` los vectores definidos en
n oM −1 `=1
(j)
(3.106) y θ` las soluciones de la ecuación (3.110).
`=1
Bibliografı́a
[1] R.P. Agarwal, Difference Equations and Inequalities. Theory, Methods and Applica-
tions, Marcel Dekker, New York, 1992.
[2] M.H. Alexander and D.E. Manolopoulos, A Stable linear reference potential algorithm
for solution of the quantum close-coupled equations in molecular scattering theory, J.
Chem. Phys. 86 (1987), 2044–2050.
[3] T.M. Apostol, Mathematical Analysis, Addison-Wesley, Reading, MA, 1977.
[4] O. Axelsson, Iterative Solution Methods, Cambridge Univ. Press., Cambridge, 1996.
[5] F.L. Bauer and A.S. Householder, Absolute norms and characteristic roots, Numerische
Math. 3 (1961), 241–246 .
[6] S.L. Campbell and C.D. Meyer, Jr., Generalized Inverses of Linear Transformations,
Pitman, London, 1971.
[7] M.C. Casabán and L. Jódar, J.A. Sánchez Cano, Stable numerical solution of strongly
coupled mixed diffusion problems, Applied Math. Lett. 15, (2002), 115–119.
[8] H.H. Chiu, Theory of irreducible linear systems, Quant. Appl. Math. XXVII, (1969),
87–104.
[9] R. Courant and D. Hilbert, Methods of Mathematical Physics, Vol. II, Interscience,
New York, 1962.
[10] P.K. Das, Optical Signal Processing, Springer, New York, 1991.
[11] N. Dunford and J. Schwartz, Linear Operators, Part I, Interscience, New York, 1957.
[12] R. Godement, Course d’Algebre, Hermann, Paris, 1957.
[13] G. Golub and C.F. Van Loan, Matrix Computations, Johns Hopkins Univ. Press, Bal-
timore, M.A., 1989.
[14] K. Gopolsamy, An arms race with deteriorating armaments, Math. Biosci, 37 (1977),
191–203.
[15] A.F. Harvey, Microwave Engineering, Academic Press, New York, 1963.
[16] D.J. Higham and N.J. Higham, MATLAB guide, Society for Industrial and Applied
Mathematics (SIAM), cop. 2000, Philadelphia.
109
110 Bibliografı́a
[17] A.L. Hodgkin and A.F. Huxley, A quantitative description of membrane current and
its applications to conduction in the giant axon of Loligo, J. Physiol., Vol. 117, pp.
500–544, (1952).
[18] L. Jódar, Algebraic and differential operator equations, Linear Algebra Appl., 102
(1988), 35–53.
[19] L. Jódar and M.C. Casabán, Convergent matrix pencils, Applied Math. Lett. 14 (5),
549–551, (2001).
[20] L. Jódar and M.C. Casabán, Convergent discrete numerical solution of coupled mixed
partial differential systems, Math. Comput. Model. 34, 283–297, (2001).
[21] L. Jódar and M.C. Casabán, Convergent discrete numerical solutions of strongly
coupled mixed parabolic systems, Utilitas Math. Paper (2002), (por aparecer).
[22] L. Jódar, E. Navarro On complete sets of solutions of polynomial matrix equations,

Applied Math. Lett. Vol. 3, No. 1 (1990), 15–18.
[23] L. Jódar, E. Navarro and J.A. Martı́n, Exact and analytic-numerical solutions of stron-
gly coupled mixed diffusion problems, Proc. Edinburgh Math. Soc. 43 (2000), 269–293.
[24] L. Jódar and J.A. Sánchez, Constructive stable numerical solutions of strongly coupled
diffusion mixed partial differential systems, Int. J. Computer Math. 73 No.2, 225–242,
(1999).
[25] L. Jódar, J.A. Sánchez and M.V. Ferrer, Stable numerical solutions of strongly coupled
mixed diffusion problems, Computers Math. Appl. 39 No. 1–2 (2000), 169–182.
[26] L. Jódar and J.A. Sánchez Cano, Discrete numerical solution of strongly coupled mixed
diffusion problems, Computers Math. Appl. 40 (2000), 471–489.
[27] R.J. Kee and L.R. Petzold, A differential/algebraic equation formulation of the method
of lines solution to systems of partial differential equations, Tech. Rep. SAND86–8893,
Sandia National Laboratories, 1986.
[28] M.P. Legua Fernández, Ecuaciones matriciales diferenciales, en diferencias y en deri-

vadas parciales, Tesis Doctoral, Universidad Politécnica de Valencia, 1990.
[29] J.A. Martı́n, E. Navarro, L. Jódar, Exact and analytic-numerical solutions of strongly
coupled mixed diffusion problems, Proc. Edinburg Math. Soc. 43 (2000), 1–25.
[30] M. Mascagni, A critical-boundary value problem of physiological significance for equa-

tions of nerve conduction, Com. Pure Appl. Math., XLII (1989), 213–227.
[31] V.S. Melezhik, I.V. Puzhynin, T.P. Puzhynina and L.N. Somov, Numerical solution of
a system of integrodifferential equations arising from the quantum-mechanical three-
body problem with Coulomb interaction, J. Comput. Phys. 54, 221–236, (1984).
[32] A.C. Metaxas, R.J. Meredith, Industrial Microwave Heating, Peter Peregrinus, Lon-
don, 1983.
Bibliografı́a 111
[33] L. Mirsky, An Introduction to Linear Algebra, Dover, New York, 1990.
[34] J.K. Mitchell, Conduction phenomena: from theory to geotechnical practice, Geotech-
nique, Vol. 41, No. 3 (1991), 229-340.
[35] H. Morimoto, Stability in the wave equation coupled with heat flows, Numer. Math. 4
(1962), 136–145.
[36] P.M. Morse and H. Feshbach, Methods of Theoretical Physics, McGraw-Hill, Tokio,
1953.
[37] E. Navarro, E. Ponsoda and L. Jódar, A matrix approach to the analytic-numerical

solution of mixed partial differential systems, Computers Math. Appl. 30 (1), 99–109,
(1995).
[38] J.M. Ortega, Numerical Analysis. A Second Course. Academic Press, New York, 1972.
[39] C.R. Rao and S.K. Mitra, Generalized Inverse of Matrices and its Applications, John-
Wiley, New York, 1971.
[40] S. Saks, A. Zygmund, Analytic Functions, Elsevier, Amsterdam, 1971.
[41] J.A. Sánchez Cano, Construcción de soluciones numéricas estables de problemas de

difusión fuertemente acoplados, Tesis Doctoral, Universidad Politécnica de Valencia,
2000.
[42] D. Sheng and K. Axelsson, Uncoupling of coupled flow in soil. A finite element method,
Int. J. Numer. and Analytical Meth. in Geomechanics, Vol. 19, pp. 537–553, (1995).
[43] G.D. Smith, Numerical Solution of Partial Differential Equations: Difference Methods,
Claredon, Oxford, 1978.
[44] G.W. Stewart and Ji-Guang Sun, Matrix Perturbation Theory, Academic Press, New
York, 1990.
[45] J.C. Strikwerda, Finite Difference Schemes and Partial Differential Equations, Chap-
man and Hall, New York, 1990.
[46] J.R. Weston and M. Stacey, Fusion Plasma Analysis, John Wiley ans Sons, New York,
1981.
[47] A.T. Winfree, When Times Breaks Down, Princeton Univ. Princeton, 1987.
[48] A.T. Winfree, Heart Muscle as a Reaction-Diffusion Medium: the Roles of Electrical
Potential Diffusion, Activation Front Curvature, and Anisotropy, Int. J. Bif. Chaos 7,
487–526, (March 1997).
[49] A.T. Yeung, J.K. Mitchell, Coupled fluid, electrical and chemical flows in soil, Geoth-
chnique vol. 43, No. 1 (1993), 121–134.

Sistemas Acoplados Edp PDF

Cargado por

Copyright:

Formatos disponibles

Sistemas Acoplados Edp PDF

Cargado por

Información del documento

Título original

Derechos de autor

Formatos disponibles

Compartir este documento

Compartir o incrustar documentos

Opciones para compartir

¿Le pareció útil este documento?

¿Este contenido es inapropiado?

Copyright:

Formatos disponibles

Sistemas Acoplados Edp PDF

Cargado por

Copyright:

Formatos disponibles

UNIVERSIDAD POLITÉCNICA DE VALENCIA

DEPARTAMENTO DE MATEMÁTICA APLICADA

Valencia, abril de 2002

Valencia, abril de 2002

Fdo: Lucas Jódar Sánchez

“ No existe problema alguno en ingenierı́a

Esta memoria trata sobre la construcción de soluciones numéricas estables de sistemas

Mediante la aplicación de un método de separación de variables discreto, la solución numéri-

Nuestro enfoque metodológico es alternativo frente al tratamiento algebraico más tradicio-

Aquesta memòria tracta sobre la construcció de solucions numèriques estables de sistemes

Mitjançant l’aplicació d’un mètode de separació de variables discret, la solució numèrica

By the application of a discrete separation of variables method, the proposed numerical

Our methodological approach is alternative to the more traditional algebraic treatment

1. Introducción, motivación y preliminares 1

1.1. Introducción y motivación

Debido a la complejidad que involucra el acoplamiento de los sistemas matriciales, los

La memoria está estructurada en dos capı́tulos. En el capı́tulo 2 se construyen soluciones

siendo A, B, C, D, Ak , Bk , k = 1, 2, matrices cuadradas en Cs×s , y la incógnita u y F

Las soluciones numéricas discretas para la ecuación

ut (x, t) − Auxx (x, t) = 0 , (1.3)

siendo A, B, C matrices cuadradas en Cs×s , y la incógnita u y F , V vectores en Cs .

Definición 1.2 Si v es un vector en Cm , su norma infinito la denotamos por kvk∞ , y

máx {|z| ; z ∈ σ(A)} ,

lo denotaremos por ρ(A).

En este trabajo utilizaremos un tipo de inversa generalizada, la inversa generalizada de

Definición 1.5 Si A ∈ Cm×n , entonces su inversa generalizada Moore-Penrose se denota

donde AH denota la traspuesta conjugada de A.

es compatible (i.e. tiene solución) si y sólo si

Además en este caso, la solución general de Ax = b viene dada por

donde z es un vector arbitrario de Cn .

Definición 1.6 [13, pág. 311]. Un subespacio E de Cm es invariante por la matriz A ∈

Observación 1.2 Si tenemos v ∈ Cm tal que v ∈ Nuc B y además Nuc B es un subespacio

LEMA 1.2 [23]. Sean M y N matrices en Cs×s , entonces:

i.e., d ∈ Nuc N I − M † M . Aplicando nuevamente el teorema 2.3.2 obtenemos:

i.e.,   h  i† h  i

De aquı́ multiplicando por M y N respectivamente y utilizando

por ser M † la inversa de Moore-Penrose, se sigue:

Definición 1.7 Toda matriz compleja A en Cm×m se puede descomponer de la siguiente

TEOREMA 1.2 [4]. Para una matriz A arbitraria se verifica

COROLARIO 1.1 [4]. Si P es una matriz hermı́tica (o real, simétrica) entonces

Observación 1.4 [4]. Para una matriz arbitraria se verifica

TEOREMA 1.4 [44, pág. 25]. Sean A y B matrices hermı́ticas, n × n, y A e = A + B.

Si A es una matriz hermı́tica denotaremos

λmı́n (A) = mı́n (σ(A)) y λmáx (A) = máx (σ(A)) .

y en la unión de los discos

TEOREMA 1.6 (Teorema de la Aplicación Espectral) [11].

LEMA 1.3 (Lema de Banach) [13]. Si P y Q son matrices de Cn×n y P es invertible,

Q es invertible , y P −1 − Q−1 < kP kkQkkP − Qk .

donde la submatriz E es invertible y H cuadrada ( E y H pueden ser de tamaños distintos),

Ya que E es invertible, se tiene que

Ahora trasladamos ésto a un sistema de ecuaciones en derivadas parciales de 2o orden con

Definición 1.9 Un polinomio p(λ), no nulo, se dice que es un polinomio anulador de A

TEOREMA 1.7 (Teorema de Cayley-Hamilton) [33, pág. 206].

TEOREMA 1.9 (Teorema de descomposición) [12, pág. 538].

COROLARIO 1.2 [40,pág. 76]. En particular, si P es invertible, entonces σ(P ) está en

ESQUEMAS EN DIFERENCIAS FINITAS [45, pág. 13].

Definición 1.12 (Esquema en diferencias convergente) Sea P u = f una ecuación

Pero como probar la convergencia de un esquema en general no es fácil, se utilizan dos

Definición 1.13 (Esquema en diferencias consistente)

i.e., h i† h i

COROLARIO 1.2 [40,pág. 76]. En particular, si P es invertible, entonces σ(P ) está en