Edit: Moved here from this topic: On 2D Lagrange polynomials
Related to this topic, I am trying to do interpolation to a field. I called legendre polynomial function in poly.py to construct my polynomial space of each element. It is okay for the most elements but to some elements, the condition number is still very big which leads to very bad interpolation results. As illustrated in this picture: The regions which blow up the color map are elements which have big condition numbers.
To reach this step, I interpolated mesh points to GL point and use ortho_basis_at this function to construct polynomial space. Could you point out what could possibly the reason of this? Could this happen even with orthogonal basis? And also I wonder what jac_ortho_basis_at_py this function actually calculated? I am a bit confused since it sets alpha and beta to zero like legendre polynomials.
I can’t really follow what you are doing and what you’re trying to acheive. Can you give me a bit more information?
Ah sorry. I am trying to do interpolation to the field with a new set of mesh. As my understanding, one can construct a link between solution and mesh using this formulation:
solution = Polynomial-space-at-element @ coefficient-matrix
To have this coefficient matrix, one can solve or simply invert this polynomial-space matrix. Once one can find a new mesh point is located in one background mesh element, one can construct a new polynomial-space at this point and multiply pre-calculated coefficient matrix to have solution at this point.
That’s what I did in the picture showed at the last post. I used ortho-basis-at-py to construct legendre-polynomial space. This works for the most elements, but for some elements which polynomial space matrix has very big conditional number (and if using numpy matrix rank which uses svd method to calculate rank, it will show these matrices are not full rank), solving these matrices leads to least square solution which could blow up the result directly. Does these operations make sense to you? If so, why orthogonal basis can lead to matrix not full rank. Another question while programming is that I have seen a function called jac_ortho_basis_at (PyFR/polys.py at a2d0be8fa7e20e9fe94cd0d007996e1c3a2bbfc8 · PyFR/PyFR · GitHub). I wonder what does this function calculated?
Hope this explains better.
For the correct approach please see the point sampler plugin, although note that this code is not robust enough for full grid interpolation.
As my understanding, this plugin extract the closest solution points to sample points and we do linear interpolation afterwards right? Is there anyway to improve accuracy since we are on spectral method already.
The plugin used to extract the nearest solution point. Now, it does full interpolation which is not trivial as for each point you want to query in physical space you need to first find the associated element and then find the place to query that element in reference space. This mapping is non-linear (at least for curved elements) and so Newton iterations are required.
Ah thanks. I was looking at the old version of PyFR. I will switch to a newer version.