error with some demos ( + a question on diagnostics)

[francispoulin]

Hello Mikael,

I am looking at a couple of examples to figure out how to do diagnostics better than what I'm doing and I noticed that the two Ginzburg_Laudau examples don't work on my latest install. Also, when I tried biharmonic1D.py.

Below you will see that I did update shenfun and then you will see the three sets of errors that I have.

I will point out tha tthe tests all pass so I thought my system was good to go but given these errors maybe there is a problem.

After I can get these working I am curious to ask your advice on how best to have diagnostics in shenfun.

Cheers!

P.S. I see that I didn't have mayavi on my system so I tried to install it using pip3 but it failed to build. The essence seems to be in this statement here. Sadly, I still get an error tha tmayavi_show is not defined.

Failed to build mayavi
Installing collected packages: configobj, traits, pyface, traitsui, apptools, envisage, pygments, vtk, mayavi
    Running setup.py install for mayavi ... done
  DEPRECATION: mayavi was installed using the legacy 'setup.py install' method, because a wheel could not be built for it. pip 21.0 will remove support for this functionality. A possible replacement is to fix the wheel build issue reported above. You can find discussion regarding this at https://github.com/pypa/pip/issues/8368.

Downloading and Extracting Packages
pip-20.2.3           | 1.7 MB    | ############################################################################### | 100% 
setuptools-49.6.0    | 742 KB    | ############################################################################### | 100% 
attrs-20.2.0         | 42 KB     | ############################################################################### | 100% 
openssl-1.1.1h       | 2.5 MB    | ############################################################################### | 100% 
Preparing transaction: done
Verifying transaction: done
Executing transaction: done

(shenfun) fpoulin@vortex2:~/software/shenfun/demo$ python Ginzburg_Landau.py 
No wisdom imported
Traceback (most recent call last):
  File "Ginzburg_Landau.py", line 91, in <module>
    U_hat = integrator.solve(U, U_hat, dt, (0, end_time))
  File "/home/fpoulin/software/anaconda3/envs/shenfun/lib/python3.8/site-packages/shenfun/utilities/integrators.py", line 241, in solve
    self.U_hat0[:] = u_hat*self.ehL_h
ValueError: operands could not be broadcast together with shapes (129,129) (2,) 

(shenfun) fpoulin@vortex2:~/software/shenfun/demo$ python Ginzburg_Landau_sphere_IRK3.py 
Traceback (most recent call last):
  File "Ginzburg_Landau_sphere_IRK3.py", line 10, in <module>
    from mayavi import mlab
ModuleNotFoundError: No module named 'mayavi'
(

(shenfun) fpoulin@vortex2:~/software/shenfun/demo$ python biharmonic1D.py 24 chebyshev
Error=6.1129316770884646e-05
Traceback (most recent call last):
  File "biharmonic1D.py", line 88, in <module>
    assert np.linalg.norm(uj-uq) < 1e-8
AssertionError

[mikaem]

There are not tests for all demos, so perhaps Ginzburg Landau is broken. Should hopefully be easy to fix😀

Try biharmonic with 36 instead of 24😀

I think I installed mayavi from conda-forge.

[francispoulin]

I suspect that Ginzburg Landau isn't too broken as it worked a few months ago and I suspect it's not a big deal.

Thanks for the suggestion on biharmonic and see that it does work. The problem was not with the code running, but the error wasn't small enough. Clearly I was being too quick. Sorry.

I uninstalled mayavi from pip and try and instead used conda-forge. It did not complain but gave Generic Warnings and then a seg fault. This isn't a shenfun library so I suppose I should contact the developers to see why it doesn't work.

[francispoulin]

Question on diagnostics:

In my simulations I want to measure things like energy at each time step, and a few other scalar fields.

I am trying to write those into a file using HDF5, which is maybe overkill but thought it should work. Unfortuantely, it does no seem to like to write scalars.

What would you recommend I use in that library to write a scalar?

[mikaem]

I think I would store the energy as a 1D dataset and keep the energy in python in a 1D numpy array.

[francispoulin]

Thanks for the advice. I'm trying to do this but having some difficulties with the syntax.

Below I do 2 things. First I write a 2D array and that works out very well. Second, I have a 1D array that has 2 elements, energy plus two other scalarys, and I want to write that into a different file and that fails, with an error copied below.

I've tried a bunch of different things but have not been able to figure it out. Might you be able to tell me the correct syntax?

    file = HDF5File("QG_1Layer_{}.h5".format(N[0]),
                     mode='w',
                     domain=((0,1),(0,1)) )
    file.write(tstep,  {'q': [ q ] } )

    norms = np.zeros(3)
    file_d = HDF5File("QG_1Layer_{}_diag.h5".format(N[0]),
                     mode='w',
                     domain=((0,1)) )
    file_d.write(tstep, {'norms': [ norms ]} )

Traceback (most recent call last):
  File "qg_1L_pert_shenfun.py", line 191, in <module>
    q_hat = integrator.solve(q, q_hat, dt, t_domain)
  File "/home/fpoulin/software/anaconda3/envs/shenfun/lib/python3.8/site-packages/shenfun/utilities/integrators.py", line 321, in solve
    self.update(u, u_hat, t, tstep, **self.params)
  File "qg_1L_pert_shenfun.py", line 152, in update
    file_d.write(tstep, {'norms': [ norms ]} )
  File "/home/fpoulin/software/anaconda3/envs/shenfun/lib/python3.8/site-packages/mpi4py_fft/io/h5py_file.py", line 118, in write
    FileBase.write(self, step, fields, **kw)
  File "/home/fpoulin/software/anaconda3/envs/shenfun/lib/python3.8/site-packages/mpi4py_fft/io/file_base.py", line 65, in write
    self._check_domain(group, u)
  File "/home/fpoulin/software/anaconda3/envs/shenfun/lib/python3.8/site-packages/mpi4py_fft/io/h5py_file.py", line 39, in _check_domain
    assert len(self.domain) == field.dimensions
AttributeError: 'numpy.ndarray' object has no attribute 'dimensions'

[mikaem]

You cannot use the write method to do this. That is for storing distributed arrays with a mesh. Open the hdf5 file, create dataset and store your 1d array as I'm sure you can find out how in the h5py documentation😀

[francispoulin]

Ah, so I've basically been trying to fit a square peg in a round hole. I guess I should be glad it didn't work out in the end.

Thanks for clarifying and will use a simple hdf5 file for the diagnostics, as that will be possible.

[francispoulin]

this is closed but I have a question on diagnostics so I thought why not post it here?

I am basing my code on Ginzburg_Landau.py and having a problem with the update function.

In the update function I am writing diagnostics, energy and the like. The variable u is part of the arguments so I thought it would be have the current value but it only seems to update it every plot or write step. It still writes the energy but it's the same energy at every plot or write step. Is this expected to happen?

Unfortuantely, this function doesn't run yet so I can't modify it easily to demosntrate what I'm trying to say. Sorry.

[mikaem]

The variable u is part of the arguments so I thought it would be have the current value but it only seems to update it every plot or write step.

U is not part of the solver, because we solve in spectral space for U_hat and not U. So either compute the energy from U_hat or transform backward to U, like in every write step or plot step.

[francispoulin]

Ah, yes, that would explain it. Thanks for explaining that and it now makes sense.