Rory/ab13md

[roryyorke]

I think this is fine as-is; some remarks (or, depending on your point of view, gripes):

1. Do we need to pass n ? Why don't we get it from Z.shape[0] ? I followed the convention of other Slycot routines, but I don't see why we do this.
2. Is there a reason we don't call np.atleast_1d and np.atleast_2d as necessary on Slycot function args? This would improve caller-friendliness, and, IMO, match expectations for Numpy-associated code.
3. The workspace arrays are optimally sized using results from Lapack's ILAENV, but I don't see that exposed in scipy.linalg.lapack. I didn't investigate the option of calling AB13MD to query optimal size parameters. I don't know if this is a major issue.
4. On the scalar real case: I've opened Structured singular value of real scalar with real uncertainty SLICOT/SLICOT-Reference#4; assuming this is actually a bug, this case could be detected in Slycot, but I'm inclined to rather handle it in python-control if I get around to adding a wrapper there.

[roryyorke]

I don't understand the MacOS failures.

From one of the the failing runs, there's something that looks weird, though it may be unrelated: https://github.com/python-control/Slycot/runs/6246110913?check_suite_focus=true#step:4:1970

Processing dependencies for slycot==0.4.0
Searching for slycot==0.4.0
Reading https://pypi.org/simple/slycot/
Downloading https://files.pythonhosted.org/packages/85/21/4e7110462f3529b2fbcff8a519b61bf64e0604b8fcbe9a07649c9bed9d7a/slycot-0.4.0.0.tar.gz#sha256=1d9921e9b04a5b9892870fd3481f7b08e6fa083a1a3848ad262819de19eb5e02

I assume this is a mistake -- we don't want to download slycot 0.4.0 from PyPI for testing, do we?

I don't follow the final error message; I guess it failed, but why?

Running slycot-0.4.0.0/setup.py -q bdist_egg --dist-dir /tmp/easy_install-w0mc9zxi/slycot-0.4.0.0/egg-dist-tmp-qcfmz0w8
    sys.exit(main())
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/cli/main_build.py", line 488, in main
    execute(sys.argv[1:])
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/cli/main_build.py", line 477, in execute
    outputs = api.build(args.recipe, post=args.post, test_run_post=args.test_run_post,
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/api.py", line 186, in build
    return build_tree(
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/build.py", line 3088, in build_tree
    packages_from_this = build(metadata, stats,
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/build.py", line 2211, in build
    utils.check_call_env(cmd, env=env, rewrite_stdout_env=rewrite_env,
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/utils.py", line 410, in check_call_env
    return _func_defaulting_env_to_os_environ('call', *popenargs, **kwargs)
  File "/usr/local/miniconda/envs/build-env/lib/python3.9/site-packages/conda_build/utils.py", line 390, in _func_defaulting_env_to_os_environ
    raise subprocess.CalledProcessError(proc.returncode, _args)
subprocess.CalledProcessError: Command '['/bin/bash', '-o', 'errexit', '/usr/local/miniconda/envs/build-env/conda-bld/slycot_1651398745184/work/conda_build.sh']' returned non-zero exit status 1.
Error: Process completed with exit code 1.

[bnavigator]

1. Do we need to pass n ? Why don't we get it from Z.shape[0] ? I followed the convention of other Slycot routines, but I
   don't see why we do this.

Some do, some don't. See e.g. mb05md and mb05nd. Changing API is a hassle but when introducing new wrappers, I would also prefer inferring n from the array shapes.

2. Is there a reason we don't call np.atleast_1d and np.atleast_2d as necessary on Slycot function args? This would improve caller-friendliness, and, IMO, match expectations for Numpy-associated code.

Not sure if this would by better to apply in one wrapper layer further up (python-control). The docstring says, give by an array, so it expects an array.

3. The workspace arrays are optimally sized using results from Lapack's ILAENV, but I don't see that exposed in scipy.linalg.lapack. I didn't investigate the option of calling AB13MD to query optimal size parameters. I don't know if this is a major issue.

If I understand correctly, it only makes a difference in the optimal size for performance.

• We could wrap ILAENV it ourselves.
• The optimal workspace size will also be returned in LDWORK(1) and LZWORK(1). Could use that after the first l iteration a performance sensitive loop?

4. On the scalar real case: I've opened Structured singular value of real scalar with real uncertainty SLICOT/SLICOT-Reference#4; assuming this is actually a bug, this case could be detected in Slycot, but I'm inclined to rather handle it in python-control if I get around to adding a wrapper there.

Agreed.

[bnavigator]

fact and n are not hidden

[roryyorke]

So "hidden" here is not arguments hidden from the caller of slycot.ab13md?

[bnavigator]

Those are the arguments of the Fortran subroutine hidden in the Python call signature of _wrapper.ab13md(). The arg_list is used for the error message handler which references the argument count and order of the Fortran routine.

[roryyorke]

OK, will fix.

[roryyorke]

Hm, what about m and info? m is apparently optional, despite my attempts to make it required in analysis.pyf, and info is hidden in ab13ed.

[bnavigator]

Not sure about m. Maybe f2py is confused by the dimension(m) and depend(m) stuff.

info shouldn't be hidden. It is in the output and you need it for the error handler.

[bnavigator]

I assume this is a mistake -- we don't want to download slycot 0.4.0 from PyPI for testing, do we?

It is building 0.4.0(.0) something and the conda can't resolve the requirement. I assume some conda or scikit-build related regression. Definitely not relevant to this PR

[roryyorke]

when introducing new wrappers, I would also prefer inferring n from the array shapes.

Great, in that case I'll remove n.

Not sure if this would by better to apply in one wrapper layer further up (python-control). The docstring says, give by an array, so it expects an array.

Fair enough. Also, I feel like if I do it here I'll be expected to change the whole of Slycot to match!

If I understand correctly, it only makes a difference in the optimal size for performance.

Yes, pretty sure that's what it's for.

The optimal workspace size will also be returned in LDWORK(1) and LZWORK(1). Could use that after the first l iteration a performance sensitive loop?

That would probably be easiest. I don't think we can easily wrap ILAENV -- it's in LAPACK, not SLICOT. I guess we could access via ctypes?

I don't want to do that as part of this PR, though; I'll leave it for when someone is working with big enough systems that it makes a difference.

I think this could be added in a backward-compatible way: add keyword-only arguments for optimal parameters as a tuple, and add a separate Python function (get_ab13md_opt) to query them. The .pyf interface would have to change to expose ldwork etc.

[bnavigator]

I don't think we can easily wrap ILAENV -- it's in LAPACK, not SLICOT. I guess we could access via ctypes?

You can provide f2py with a Fortran wrapper routine calling ILAENV. See e.g. the XERBLA override.

[roryyorke]

Thanks @bnavigator .

[bnavigator]

Unless you want to also hide m the same way as n, this LGTM.

[roryyorke]

Good idea.