Acoustic kinetic/potential energy computation

Usage
1

I discovered that the acoustic kinetic and potential energy are somewhat “interestingly” computed.

As it is not stated on the docu page how they are derived, I took a look into the source code.

  1. Acoustic kinetic energy:

    <regionResult type="acouKinEnergy"/>

I’d expect this to be something like rho0 v^2/2, which is also stated in a comment in the code. But the code says something different:

shared_ptr<BaseFeFunction> deriv1vFct;
    if ( formulation_ == ACOU_PRESSURE )
      deriv1vFct = timeDerivFeFunctions_[ACOU_PRESSURE_DERIV_1];
    else
      deriv1vFct = timeDerivFeFunctions_[ACOU_POTENTIAL_DERIV_1];

    if( isComplex_ ) {
      keFunc.reset(new EnergyResultFunctor<Complex>(deriv1vFct, kinEnergy, 0.5));
    } else {
      keFunc.reset(new EnergyResultFunctor<Double>(deriv1vFct, kinEnergy, 0.5));
    }
    resultFunctors_[ACOU_KIN_ENERGY] = keFunc;
    massFormFunctors_.insert(keFunc);

  1. Acoustic potential energy:

    <regionResult type="acouPotEnergy"/>

Here, a comment states that the quantity should be something like p^2/(2 rho0 c^2), but the code says

    if( isComplex_ ) {
      keFuncPot.reset(new EnergyResultFunctor<Complex>(feFct, potEnergy, 0.5));
    } else {
      keFuncPot.reset(new EnergyResultFunctor<Double>(feFct, potEnergy, 0.5));
    }
    resultFunctors_[ACOU_POT_ENERGY] = keFuncPot;
    stiffFormFunctors_.insert(keFuncPot);

Interpretation

The kinetic energy does not make sense to me and appears like it has been taken from the mechPDE, where the time derivative is required to obtain the velocity. However, in the acousticPDE, we would need something with a gradient operation.

The potential energy appears somewhat correct for the pressure formulation, as the feFct is directly treated. For the potential formulation, however, a time derivative would be required.

For both, I am not sure where and if the weighting by material parameters applies. Possibly this is included in the EnergyResultFunctor, but I didn’t find time to look more deeply into its implementation.

@ftoth @dmayrhofer , I guess you have ideas for this. Should I open an issue?
I guess the fix could be of low effort.

2

To be honest I have not used this result so I have never checked if it is correct. But simply from taking a look at the basically equal treatment for pressure and potential formulation, I suspect that there is a bug. So you initial idea of opening an issue is the right call imo. Do you need the result urgently? Because this could be a nice implementation issue for the Multiphysics course.

3

would be fine for my side. I’ll open the issue.

4

would be fine for my side.

I just discovered that there is already an open issue discussing the topic.

I’ll provide a link to this chat there.