Hello everyone,
I want to get the information of the adjacont element, such as rho, u… I know it’s not a simple operation, because pyfr is unstructured solver. But it has been bothering me a long time. In my opinipn, I use Mako to get the information of the adjacent element, then I can do some operation in element.py to make use of these information. Could someone give me some more details?
Second, how to get the length of the current element?
With best regards
I would suggest writing down what you want to do mathematically before trying to think about how to implement it. Lets assume you have access to the full solution polynomials of all neighbouring elements, what would you do with it? And how would these operations maintain high order accuracy? How would you account for different polynomial spaces in adjacent elements? (Tetrahedron with neighbours that are a prism, another tetrahedron, and two pyramids, for example.) All of this becomes even more complex when you remember that these polynomials are by default defined nodally and so you need to be willing to deal with different solution point types.
Assuming that everything above has been worked out the first thing you’ll need to do is handle the MPI transfers between ranks (as you’ll now need to send entire elements to neighbouring ranks).
Regards, Freddie.
At present, I only consider the rectangular element cases, it could be more easy to implement it. and how to have access to the solution polynomials of neighbouring element (the density is what I only need). Of course the aim of this operation is to implement KXRCF indicator.
Thaks a lot.
Of course I can transfer entropy instead of solution polynomials.
I have the KXRCF sensor implemented here: PyFR/pyfr/solvers/euler/kernels/kxrcf.mako at feature/cmp · WillTrojak/PyFR · GitHub
We used it in a multi-species solver work: [2308.02426v1] Positivity-preserving discontinuous spectral element methods for compressible multi-species flows
Although we only used the sensor to switch the entropy part of the entropy filter
Thanks, I implement KXRCF, and it comes the following error, could you tell me how to solve it?
Traceback (most recent call last):
File "/home/bin/pyfr", line 33, in <module>
sys.exit(load_entry_point('pyfr', 'console_scripts', 'pyfr')())
File "_main__.py", line 126, in main
args.process(args)
File "__main__.py", line 273, in process_run
_process_common(
File "__main__.py", line 269, in _process_common
solver.run()
File "integrators/base.py", line 141, in run
self.advance_to(t)
File "integrators/std/controllers.py", line 74, in advance_to
idxcurr = self.step(self.tcurr, dt)
File "integrators/std/steppers.py", line 57, in step
limit(t, r0)
File "pyfr/solvers/base/system.py", line 279, in limit
for graph in self._limit_graphs(uinbank):
File "pyfr/util.py", line 39, in newmeth
res = cache[key] = meth(self, *args, **kwargs)
File "pyfr/solvers/baseadvec/system.py", line 143, in _limit_graphs
g1.add_all(k['iint/comm_jump'], deps=k['mpiint/jump_fpts_pack'] + k['eles/disu'])
File "pyfr/backends/base/types.py", line 342, in add_all
self.add(k, deps)
File "pyfr/backends/base/types.py", line 331, in add
rdeps = [self.knodes[d] for d in deps]
File "pyfr/backends/base/types.py", line 331, in <listcomp>
rdeps = [self.knodes[d] for d in deps]
KeyError: <pyfr.backends.openmp.provider.OpenMPKernel object at 0x7fd9cf8e96c0>
This suggests you are passing an invalid kernel as a dependency (as in there is a problem with your deps argument with one or more of the kernels passed not being present in the graph).
Regards, Freddie.
Thanks for your replay.
I don’t really understand that part of the code in baseadvec/system. How to choose the kernels to pass to deps? Are there any rules for naming parameters,such as ‘eles/disu’.
And in kxrcf, have you determined which is the inflow boundary? Last but not least, I want to obtain kxrcf value in all interfaces of an element, could you give me some tips？
Best regards
You determine the dependencies based on what inputs your kernel requires. If your kernel makes use of matrix A then any kernels which write to matrix A need to be dependencies so as to ensure that they are run before your kernel.
The naming convention is mostly arbitrary. The eles/ prefix says that the kernel comes from an Elements class and disu means “discontinuous solution” since the kernel computes the discontinuous solution at the flux points.
Kernels that need to run at all interface points typically need to come in triplets. One kernel for the internal interfaces, one for boundary interfaces, and one for MPI interfaces.
Regards, Freddie.
Thanks a lot. I solved the dependency problem. But why did this error occur?
pyfr", line 33, in <module>
sys.exit(load_entry_point('pyfr', 'console_scripts', 'pyfr')())
pyfr/__main__.py", line 126, in main
args.process(args)
pyfr/__main__.py", line 273, in process_run
_process_common(
pyfr/__main__.py", line 269, in _process_common
solver.run()
pyfr/integrators/base.py", line 141, in run
self.advance_to(t)
pyfr/integrators/std/controllers.py", line 77, in advance_to
self._accept_step(dt, idxcurr)
pyfr/integrators/std/controllers.py", line 41, in _accept_step
self._run_plugins()
pyfr/integrators/base.py", line 96, in _run_plugins
plugin(self)
pyfr/plugins/residual.py", line 50, in __call__
resid = sum(norm(dt_s)**self.lp for dt_s in intg.dt_soln)
pyfr/plugins/residual.py", line 50, in <genexpr>
resid = sum(norm(dt_s)**self.lp for dt_s in intg.dt_soln)
pyfr/plugins/residual.py", line 46, in <lambda>
norm = lambda x: np.linalg.norm(x, axis=(0, 2), ord=self.lp)
pyfr/lib/python3.10/site-packages/numpy/linalg/linalg.py", line 2601, in norm
ret = _multi_svd_norm(x, row_axis, col_axis, amax)
pyfr/lib/python3.10/site-packages/numpy/linalg/linalg.py", line 2375, in _multi_svd_norm
result = op(svd(y, compute_uv=False), axis=-1)
pyfr/lib/python3.10/site-packages/numpy/linalg/linalg.py", line 1693, in svd
s = gufunc(a, signature=signature, extobj=extobj)
pyfr/lib/python3.10/site-packages/numpy/linalg/linalg.py", line 121, in _raise_linalgerror_svd_nonconvergence
raise LinAlgError("SVD did not converge")
numpy.linalg.LinAlgError: SVD did not converge
This looks like a bug in the plugin. I will try and get it fixed over the next few days.
Regards, Freddie.
But when I’m not using the implemented limiter, it works without error. Are there any possible mistakes I made?
It is possible your limiter is causing highly erratic residuals. I would suggest running without the plugin and seeing if the results are reasonable.
Regards, Freddie.
When I run without the plugin, the variable value of the result is zero. And I simulated the sod lax problem.
The residual is zero or the solution is zero?
Regards, Freddie.
The solution value is zero, I turned off the residual plugin.
This seems to suggest there is an issue with your limiter.
Regards, Freddie.
When the case is runing, the solver didn’t call the function “limiter” in euler/element. Whether graph has errors? Could you help me check it? Thanks a lot.
@memoize
def _limiter_graphs(self, uinbank):
m = self._mpireqs
k, _ = self._get_kernels(uinbank, None)
def deps(dk, *names):
return self._kdeps(k, dk, *names)
g1 = self.backend.graph()
g1.add_all(k['eles/disu'])
g1.add_mpi_reqs(m['mean_fpts_recv'])
# Compute local mean within element
g1.add_all(k['eles/element_mean'])
# If shock sensor, pack and send the jump values to neighbours
g1.add_mpi_reqs(m['jump_fpts_recv'])
#g1.add_all(k['eles/kxrcf'])
g1.add_all(k['mpiint/jump_fpts_pack'], deps=k['eles/disu'])
for send, pack in zip(m['jump_fpts_send'], k['mpiint/jump_fpts_pack']):
g1.add_mpi_req(send, deps=[pack])
# Compute sensor interface jump terms at internal/boundary interfaces
g1.add_all(k['iint/comm_jump'], deps=k['mpiint/jump_fpts_pack'] + k['eles/disu'])
g1.add_all(k['bcint/comm_jump'], deps=k['eles/disu'])
# Pack and send the mean values to neighbors
g1.add_all(k['mpiint/mean_fpts_pack'], deps=k['eles/element_mean'])
for send, pack in zip(m['mean_fpts_send'], k['mpiint/mean_fpts_pack']):
g1.add_mpi_req(send, deps=[pack])
# Compute common entropy minima at internal/boundary interfaces
g1.add_all(k['iint/comm_mean'], deps=k['eles/element_mean'])
g1.add_all(k['bcint/comm_mean'],
deps=k['eles/element_mean'] + k['eles/disu'])
g1.add_all(k['eles/kxrcf'],
deps=k['iint/comm_jump'] + k['bcint/comm_jump'])
g1.add_all(k['eles/limiter'],
deps=k['iint/comm_jump'] + k['bcint/comm_jump']
+ k['iint/comm_mean'] + k['bcint/comm_mean'])
# g1.add_all(k['eles/limiter'])
g1.commit()
#if 'mpiint/comm_jump' in k:
# Compute common entropy minima at MPI interfaces
g2 = self.backend.graph()
# Compute sensor jumps at MPI interfaces
g2.add_all(k['mpiint/jump_fpts_unpack'])
for l in k['mpiint/comm_jump']:
g2.add(l, deps=deps(l, 'mpiint/jump_fpts_unpack'))
g2.add_all(k['mpiint/mean_fpts_unpack'])
for l in k['mpiint/comm_mean']:
g2.add(l, deps=deps(l, 'mpiint/mean_fpts_unpack'))
# Compute sensor values at elements
g2.add_all(k['eles/kxrcf'], deps=k['mpiint/comm_jump'])
g2.add_all(k['eles/limiter'],
deps=k['mpiint/comm_jump'] + k['mpiint/comm_mean'])
# g2.add_all(k['eles/limiter'])
g2.commit()
return g1, g2
#else:
#return g1,s
Did the other kernels in the graph get called? Is it possible that the graph is not being run at all?
Further, can you confirm that k['eles/limiter'] is not an empty list?
Regards, Freddie.
No.
The graph is running.
The limiter knernel exists, and it is called in euler/element. This isn’t an empty list?