Skip to content

Commit d2f684d

Browse files
committed
extend the parser to warnings
1 parent 1877bc3 commit d2f684d

File tree

7 files changed

+273
-223
lines changed

7 files changed

+273
-223
lines changed

slycot/exceptions.py

Lines changed: 98 additions & 70 deletions
Original file line numberDiff line numberDiff line change
@@ -20,6 +20,8 @@
2020

2121
import re
2222

23+
from warnings import warn
24+
2325

2426
class SlycotError(RuntimeError):
2527
"""Slycot exception base class"""
@@ -45,8 +47,91 @@ class SlycotArithmeticError(SlycotError, ArithmeticError):
4547
pass
4648

4749

50+
class SlycotWarning(UserWarning):
51+
"""Slycot Warning"""
52+
53+
def __init__(self, message, info):
54+
super(SlycotWarning, self).__init__(message)
55+
self.info = info
56+
57+
58+
class SlycotResultWarning(SlycotWarning):
59+
"""Slycot computation result warning
60+
61+
A Slycot routine returned a nonzero info parameter that warns about the
62+
returned results, but the results might still be usable.
63+
"""
64+
65+
pass
66+
67+
68+
def _parse_docsection(section_name, docstring, checkvars):
69+
slycot_error = None
70+
message = None
71+
docline = iter(docstring.splitlines())
72+
try:
73+
74+
info_eval = False
75+
while section_name not in next(docline):
76+
continue
77+
section_indent = next(docline).index("-")
78+
79+
for l in docline:
80+
# ignore blank lines
81+
if not l.strip():
82+
continue
83+
84+
# reached next section without match
85+
if l[section_indent] == "-":
86+
break
87+
88+
# Exception Type
89+
ematch = re.match(
90+
r'(\s*)(Slycot.*(Error|Warning)) : e', l)
91+
if ematch:
92+
error_indent = len(ematch.group(1))
93+
slycot_error = ematch.group(2)
94+
95+
# new infospec
96+
if slycot_error:
97+
imatch = re.match(
98+
r'(\s{' + str(error_indent + 1) + r',}):(.*):\s*(.*)', l)
99+
if imatch:
100+
infospec_indent = len(imatch.group(1))
101+
infospec = imatch.group(2)
102+
# Don't handle the standard case unless we have i
103+
if infospec == "e.info = -i":
104+
if 'i' not in checkvars.keys():
105+
continue
106+
infospec_ = infospec.replace(" = ", " == ")
107+
try:
108+
info_eval = eval(infospec_, checkvars)
109+
except NameError:
110+
raise RuntimeError("Unknown variable in infospec: "
111+
+ infospec)
112+
except SyntaxError:
113+
raise RuntimeError("Invalid infospec: " + infospec)
114+
if info_eval:
115+
message = imatch.group(3).strip() + '\n'
116+
mmatch = re.match(
117+
r'(\s{' + str(infospec_indent+1) + r',})(.*)',
118+
next(docline))
119+
if not mmatch:
120+
break # docstring
121+
body_indent = len(mmatch.group(1))
122+
message += mmatch.group(2) + '\n'
123+
for l in docline:
124+
if l and not l[:body_indent].isspace():
125+
break # message body
126+
message += l[body_indent:] + '\n'
127+
break # docstring
128+
except StopIteration:
129+
pass
130+
return (slycot_error, message)
131+
132+
48133
def raise_if_slycot_error(info, arg_list=None, docstring=None, checkvars={}):
49-
"""Raise exceptions if slycot info returned is non-zero.
134+
"""Raise exceptions or warnings if slycot info returned is non-zero.
50135
51136
Parameters
52137
----------
@@ -80,6 +165,9 @@ def raise_if_slycot_error(info, arg_list=None, docstring=None, checkvars={}):
80165
a generic SlycotParameterError is raised if no custom text was defined in
81166
the docstring or no docstring is provided.
82167
168+
To rase warnings, define a "Warns" section similarly formatted as "Raises"
169+
using the ``SlycotResultWarning : e`` definition name.
170+
83171
Example
84172
-------
85173
>>> def fun(info):
@@ -111,77 +199,17 @@ def raise_if_slycot_error(info, arg_list=None, docstring=None, checkvars={}):
111199
SlycotArithmeticError: 4 is between 4 and 1.2e+02!
112200
"""
113201
if docstring:
114-
slycot_error_map = {"SlycotError": SlycotError,
115-
"SlycotParameterError": SlycotParameterError,
116-
"SlycotArithmeticError": SlycotArithmeticError}
117-
118-
docline = iter(docstring.splitlines())
119-
info_eval = False
120-
try:
121-
while "Raises" not in next(docline):
122-
continue
202+
checkvars['e'] = SlycotError("", info)
123203

124-
section_indent = next(docline).index("-")
125-
126-
slycot_error = None
127-
for l in docline:
128-
print(l)
129-
# ignore blank lines
130-
if not l.strip():
131-
continue
132-
133-
# reached end of Raises section without match
134-
if not l[:section_indent].isspace():
135-
return None
136-
137-
# Exception Type
138-
ematch = re.match(
139-
r'(\s*)(Slycot(Parameter|Arithmetic)?Error) : e', l)
140-
if ematch:
141-
error_indent = len(ematch.group(1))
142-
slycot_error = ematch.group(2)
143-
144-
# new infospec
145-
if slycot_error:
146-
imatch = re.match(
147-
r'(\s{' + str(error_indent + 1) + r',}):(.*):\s*(.*)',
148-
l)
149-
if imatch:
150-
infospec_indent = len(imatch.group(1))
151-
infospec = imatch.group(2)
152-
# Don't handle the standard case unless we have i
153-
if infospec == "e.info = -i":
154-
if 'i' not in checkvars.keys():
155-
continue
156-
infospec_ = infospec.replace(" = ", " == ")
157-
checkvars['e'] = SlycotError("", info)
158-
try:
159-
info_eval = eval(infospec_, checkvars)
160-
except NameError:
161-
raise RuntimeError("Unknown variable in infospec: "
162-
+ infospec)
163-
except SyntaxError:
164-
raise RuntimeError("Invalid infospec: "
165-
+ infospec)
166-
if info_eval:
167-
message = imatch.group(3).strip() + '\n'
168-
mmatch = re.match(
169-
r'(\s{' + str(infospec_indent+1) + r',})(.*)',
170-
next(docline))
171-
if not mmatch:
172-
break # docstring
173-
body_indent = len(mmatch.group(1))
174-
message += mmatch.group(2) + '\n'
175-
for l in docline:
176-
if l and not l[:body_indent].isspace():
177-
break # message body
178-
message += l[body_indent:] + '\n'
179-
break # docstring
180-
except StopIteration:
181-
pass
182-
if info_eval and message:
204+
exception, message = _parse_docsection("Raises", docstring, checkvars)
205+
if exception and message:
183206
fmessage = '\n' + message.format(**checkvars).strip()
184-
raise slycot_error_map[slycot_error](fmessage, info)
207+
raise globals()[exception](fmessage, info)
208+
209+
warning, message = _parse_docsection("Warns", docstring, checkvars)
210+
if warning and message:
211+
fmessage = message.format(**checkvars).strip()
212+
warn(globals()[warning](fmessage, info))
185213

186214
if info < 0 and arg_list:
187215
message = ("The following argument had an illegal value: {}"

slycot/math.py

Lines changed: 69 additions & 62 deletions
Original file line numberDiff line numberDiff line change
@@ -578,63 +578,78 @@ def mb05md(a, delta, balanc='N'):
578578
579579
Matrix exponential for a real non-defective matrix
580580
581-
To compute exp(A*delta) where A is a real N-by-N non-defective
581+
To compute ``exp(A*delta)`` where `A` is a real N-by-N non-defective
582582
matrix with real or complex eigenvalues and delta is a scalar
583583
value. The routine also returns the eigenvalues and eigenvectors
584-
of A as well as (if all eigenvalues are real) the matrix product
585-
exp(Lambda*delta) times the inverse of the eigenvector matrix of
586-
A, where Lambda is the diagonal matrix of eigenvalues.
584+
of `A` as well as (if all eigenvalues are real) the matrix product
585+
``exp(Lambda*delta)`` times the inverse of the eigenvector matrix of
586+
`A`, where `Lambda` is the diagonal matrix of eigenvalues.
587587
Optionally, the routine computes a balancing transformation to
588588
improve the conditioning of the eigenvalues and eigenvectors.
589589
590-
Required arguments:
591-
A : input rank-2 array('d') with bounds (n,n)
592-
Square matrix
593-
delta : input 'd'
594-
The scalar value delta of the problem.
595-
596-
Optional arguments:
597-
balanc : input char*1
598-
Indicates how the input matrix should be diagonally scaled
599-
to improve the conditioning of its eigenvalues as follows:
600-
= 'N': Do not diagonally scale;
601-
= 'S': Diagonally scale the matrix, i.e. replace A by
602-
D*A*D**(-1), where D is a diagonal matrix chosen
603-
to make the rows and columns of A more equal in
604-
norm. Do not permute.
590+
Parameters
591+
----------
592+
A : (n, n) array_like
593+
Square matrix
594+
delta : float
595+
The scalar value delta of the problem.
596+
balanc : {'N', 'S'}, optional
597+
Indicates how the input matrix should be diagonally scaled
598+
to improve the conditioning of its eigenvalues as follows:
599+
600+
:= 'N': Do not diagonally scale;
601+
:= 'S': Diagonally scale the matrix, i.e. replace `A` by
602+
``D*A*D**(-1)``, where `D` is a diagonal matrix chosen
603+
to make the rows and columns of A more equal in
604+
norm. Do not permute.
605605
606-
Return objects:
607-
Ar : output rank-2 array('d') with bounds (n,n)
608-
Contains the solution matrix exp(A*delta)
609-
Vr : output rank-2 array('d') with bounds (n,n)
610-
Contains the eigenvector matrix for A. If the k-th
611-
eigenvalue is real the k-th column of the eigenvector
612-
matrix holds the eigenvector corresponding to the k-th
613-
eigenvalue. Otherwise, the k-th and (k+1)-th eigenvalues
614-
form a complex conjugate pair and the k-th and (k+1)-th
615-
columns of the eigenvector matrix hold the real and
616-
imaginary parts of the eigenvectors corresponding to these
617-
eigenvalues as follows. If p and q denote the k-th and
618-
(k+1)-th columns of the eigenvector matrix, respectively,
619-
then the eigenvector corresponding to the complex
620-
eigenvalue with positive (negative) imaginary value is
621-
given by
622-
p + q*j (p - q*j), where j^2 = -1.
623-
Yr : output rank-2 array('d') with bounds (n,n)
624-
contains an intermediate result for computing the matrix
625-
exponential. Specifically, exp(A*delta) is obtained as the
626-
product V*Y, where V is the matrix stored in the leading
627-
N-by-N part of the array V. If all eigenvalues of A are
628-
real, then the leading N-by-N part of this array contains
629-
the matrix product exp(Lambda*delta) times the inverse of
630-
the (right) eigenvector matrix of A, where Lambda is the
631-
diagonal matrix of eigenvalues.
632-
633-
VAL : output rank-1 array('c') with bounds (n)
634-
Contains the eigenvalues of the matrix A. The eigenvalues
635-
are unordered except that complex conjugate pairs of values
636-
appear consecutively with the eigenvalue having positive
637-
imaginary part first.
606+
Returns
607+
-------
608+
Ar : (n, n) ndarray
609+
Contains the solution matrix ``exp(A*delta)``
610+
Vr : (n, n) ndarray
611+
Contains the eigenvector matrix for `A`. If the `k`-th
612+
eigenvalue is real the `k`-th column of the eigenvector
613+
matrix holds the eigenvector corresponding to the `k`-th
614+
eigenvalue. Otherwise, the `k`-th and `(k+1)`-th eigenvalues
615+
form a complex conjugate pair and the k-th and `(k+1)`-th
616+
columns of the eigenvector matrix hold the real and
617+
imaginary parts of the eigenvectors corresponding to these
618+
eigenvalues as follows. If `p` and `q` denote the `k`-th and
619+
`(k+1)`-th columns of the eigenvector matrix, respectively,
620+
then the eigenvector corresponding to the complex
621+
eigenvalue with positive (negative) imaginary value is
622+
given by
623+
``p + q*j (p - q*j), where j^2 = -1.``
624+
Yr : (n, n) ndarray
625+
contains an intermediate result for computing the matrix
626+
exponential. Specifically, ``exp(A*delta)`` is obtained as the
627+
product ``V*Y``, where `V` is the matrix stored in the leading
628+
`n`-by-`n` part of the array `V`. If all eigenvalues of `A` are
629+
real, then the leading `n`-by-`n` part of this array contains
630+
the matrix product ``exp(Lambda*delta)`` times the inverse of
631+
the (right) eigenvector matrix of `A`, where `Lambda` is the
632+
diagonal matrix of eigenvalues.
633+
634+
VAL : (n,) real or complex ndarray
635+
Contains the eigenvalues of the matrix `A`. The eigenvalues
636+
are unordered except that complex conjugate pairs of values
637+
appear consecutively with the eigenvalue having positive
638+
imaginary part first.
639+
640+
Warns
641+
------
642+
SlycotResultWarning : e
643+
:0 < e.info <=n:
644+
the QR algorithm failed to compute all
645+
the eigenvalues; no eigenvectors have been computed;
646+
w[{e.info:n}] contains eigenvalues which have converged;
647+
:e.info == n+1:
648+
The inverse of the eigenvector matrix could not
649+
be formed due to an attempt to divide by zero, i.e.,
650+
the eigenvector matrix is singular;
651+
:e.info == n+2:
652+
Matrix A is defective, possibly due to rounding errors.
638653
"""
639654
hidden = ' (hidden by the wrapper)'
640655
arg_list = ['balanc', 'n', 'delta', 'a', 'lda'+hidden, 'v', 'ldv'+hidden,
@@ -647,17 +662,9 @@ def mb05md(a, delta, balanc='N'):
647662
delta=delta,
648663
a=a)
649664

650-
raise_if_slycot_error(INFO, arg_list)
651-
652-
if INFO > 0 and INFO <= n:
653-
raise SlycotArithmeticError("Incomplete eigenvalue calculation, "
654-
"missing {} eigenvalues".format(INFO),
655-
INFO)
656-
elif INFO == n+1:
657-
raise SlycotArithmeticError("Eigenvector matrix singular", INFO)
658-
elif INFO == n+2:
659-
raise SlycotArithmeticError("Matrix A is defective, "
660-
"possibly due to rounding errors.", INFO)
665+
raise_if_slycot_error(INFO, arg_list,
666+
docstring=mb05md.__doc__, checkvars=locals())
667+
661668
if not all(VALi == 0):
662669
VAL = VALr + 1J*VALi
663670
else:

0 commit comments

Comments
 (0)