AAAReductor #2249
Draft
AAAReductor #2249
Conversation
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
artpelling
requested changes
Jan 16, 2024
src/pymor/reductors/aaa.py
Outdated
| # Define ROM with constant output | ||
| tf = np.vectorize(lambda s: np.mean(self.Hs)) | ||
|
|
||
| self.itpl_part = [] |
There was a problem hiding this comment.
Suggested change
| self.itpl_part = [] |
This line should be deleted, right?
src/pymor/reductors/aaa.py
Outdated
|
|
||
| while len(self.itpl_part) < max_itpl: | ||
|
|
||
| # compute approximation error over LS partiation of sampled data |
There was a problem hiding this comment.
Suggested change
| # compute approximation error over LS partiation of sampled data | |
| # compute approximation error over LS partition of sampled data |
| @@ -451,3 +452,239 @@ def bary_func(*args): | |||
| return np.atleast_2d(nd) | |||
|
|
|||
| return bary_func | |||
|
|
|||
src/pymor/reductors/aaa.py
Outdated
| w /= np.linalg.norm(w) | ||
| v /= np.linalg.norm(v) | ||
| self.MIMO_Hs = Hs | ||
| self.Hs = Hs @ v @ w |
There was a problem hiding this comment.
is this correct? it looks a bit dangerous..
There was a problem hiding this comment.
Yes, it is correct. This is just the "cheap" way of dealing with MIMO data in rational approximation. Implementing the actual tangential AAA is the way to go in the future, but for now I think it works well.
src/pymor/reductors/aaa.py
Outdated
Comment on lines
504
to
519
| if len(Hs.shape) > 1: | ||
| self._dim_output = Hs.shape[1] | ||
| self._dim_input = Hs.shape[2] | ||
| rng = new_rng(0) | ||
| w = rng.normal(size=(self._dim_output,)) | ||
| v = rng.normal(size=(self._dim_input,)) | ||
| w /= np.linalg.norm(w) | ||
| v /= np.linalg.norm(v) | ||
| self.MIMO_Hs = Hs | ||
| self.Hs = Hs @ v @ w | ||
| else: | ||
| self._dim_input = 1 | ||
| self._dim_output = 1 | ||
|
|
||
| if self._dim_input == self._dim_output == 1: | ||
| Hs = np.squeeze(Hs) |
There was a problem hiding this comment.
Suggested change
| if len(Hs.shape) > 1: | |
| self._dim_output = Hs.shape[1] | |
| self._dim_input = Hs.shape[2] | |
| rng = new_rng(0) | |
| w = rng.normal(size=(self._dim_output,)) | |
| v = rng.normal(size=(self._dim_input,)) | |
| w /= np.linalg.norm(w) | |
| v /= np.linalg.norm(v) | |
| self.MIMO_Hs = Hs | |
| self.Hs = Hs @ v @ w | |
| else: | |
| self._dim_input = 1 | |
| self._dim_output = 1 | |
| if self._dim_input == self._dim_output == 1: | |
| Hs = np.squeeze(Hs) | |
| if len(Hs.shape) > 1: | |
| self._dim_output = Hs.shape[1] | |
| self._dim_input = Hs.shape[2] | |
| if self._dim_input == self._dim_output == 1: | |
| Hs = np.squeeze(Hs) | |
| else: | |
| rng = new_rng(0) | |
| w = rng.normal(size=(self._dim_output,)) | |
| v = rng.normal(size=(self._dim_input,)) | |
| w /= np.linalg.norm(w) | |
| v /= np.linalg.norm(v) | |
| self.MIMO_Hs = Hs | |
| self.Hs = Hs @ v @ w |
src/pymor/reductors/aaa.py
Outdated
Comment on lines
531
to
535
| Initial partition for interpolation values. Should be `None` or a nested list | ||
| such that `itpl_part[i]` corresponds to indices of interpolated values with | ||
| respect to the `i`-th variable. I.e., `self.sampling_values[i][itpl_part[i]]` | ||
| represents a list of all initially interpolated samples of the `i`-th variable. | ||
| If `None` p-AAA will start with no interpolated values. |
src/pymor/reductors/aaa.py
Outdated
Comment on lines
537
to
539
| Maximum number of interpolation points to use with respect to each | ||
| variable. Should be `None` or a list such that `self.num_vars == len(max_itpl)`. | ||
| If `None` `max_itpl[i]` will be set to `len(self.sampling_values[i]) - 1`. |
src/pymor/reductors/aaa.py
Outdated
Comment on lines
574
to
576
| errs = np.empty(len(ls_part)) | ||
| for i, ls_s in enumerate(self.s[ls_part]): | ||
| errs[i] = np.abs(tf(ls_s) - self.Hs[ls_part[i]]) |
There was a problem hiding this comment.
Suggested change
| errs = np.empty(len(ls_part)) | |
| for i, ls_s in enumerate(self.s[ls_part]): | |
| errs[i] = np.abs(tf(ls_s) - self.Hs[ls_part[i]]) | |
| errs = np.array([np.abs(tf(self.s[i]) - self.Hs[i]) for i in ls_part]) |
| Cr = np.empty((self._dim_output,r*self._dim_input), dtype=itpl_vals.dtype) | ||
| for i in range(r): | ||
| Cr[:,self._dim_input*i:self._dim_input*(i+1)] = itpl_vals[i,:] | ||
| L = np.diag(itpl_nodes) |
There was a problem hiding this comment.
I get an error that itpl_nodes is undefined, if it converges at the first step. Maybe we can handle this edge case.
There was a problem hiding this comment.
I decided to force the algorithm to pick at least one interpolation point. Otherwise the resulting TF is just constant, thus not strictly proper and things are not working out nicely.
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
This PR adds the
AAAReductor. Although we already have thePAAAReductorwhich can also handle non-parametric problems, the newAAAReductordoes two things differently: First, it is based on strictly proper barycentric forms unlike the p-AAA implementation which uses just proper barycentric forms. Additionally, theAAAReductorallows for directly computing the ROM as anLTIModelrather than aTransferFunction.