Potential Energy and Kinteic Energy use different unit?

Message

bzhang · #1 Unread post by **bzhang** » Thu Sep 28, 2006 6:59 pm

In diag.F

          cff=g/rho0
          DO k=N(ng),1,-1
            DO i=Istr,Iend
              ke2d(i,j)=ke2d(i,j)+                                      &
     &                  Hz(i,j,k)*(rho(i,j,k)+1000.0_r8)*               &
     &                   0.5_r8*(u(i  ,j,k,nstp)*u(i,j,k,nstp)+         &
     &                           u(i+1,j,k,nstp)*u(i+1,j,k,nstp)+       &
     &                           v(i,j  ,k,nstp)*v(i,j  ,k,nstp)+       &
     &                           v(i,j+1,k,nstp)*v(i,j+1,k,nstp))
             pe2d(i,j)=pe2d(i,j)+                                       &
     &                 cff*Hz(i,j,k)*(rho(i,j,k)+1000.0_r8)*            &
     &                 (z_r(i,j,k)-z_w(i,j,0))
            END DO
          END DO

The inconsistency is from the cff term. pe2d is divided by rho0, but ke2d is not. Thus the calculation of total energy seems wrong: avgkp=avgke+avgpe without further rho0 corrections.

shchepet · #2 Unread post by **shchepet** » Fri Sep 29, 2006 5:47 am

There is an error here: multiplier *(rho(i,j,k)+1000.0_r8) should not be present in ke2d calculation. Just change it

ke2d(i,j) = ke2d(i,j) + Hz(i,j,k*(0.5_r8*(u(i ,j,k,nstp)*u(i,j,k,nstp)+....

You may also look at diag.F routine from

http://www.atmos.ucla.edu/~alex/ROMS/roms.tar

to use it for a reference.

As a matter of fact, ROMS always been a Boussinesq approximation model (unless Hernan
changed it it ROMS 3.0, which is unlikely). Consistently with Boussinesq approximation
primitive equations, KE should use CONSTANT (Boussinesq reference) density, rather than
in-situ.

....it looks like somebody wanted to improve accuracy of KE recently without taking into consideraton the rest of the code.

shchepet · #3 Unread post by **shchepet** » Fri Sep 29, 2006 5:21 pm

...something else caught my attention. It looks like there is a missing factor of 1/2
in KE computation, since

KE = intergal of V^2/2

and not just of V^2. The factor of 0.5 in ke2d computation comes from averaging,
0.5*(u(i+1,j)^2+u(i,j)^2) approximates u^2 and not u^2/2. So, unless ke2d is renormalized
later in the code, the factor of 0.5 should be changed to 0.25.

bzhang · #4 Unread post by **bzhang** » Fri Sep 29, 2006 9:37 pm

The original code about ke2d in roms-2.2 is:

Code: Select all

DO i=Istr,Iend
            ke2d(i,j)=(zeta(i,j,krhs)+h(i,j))*                          &
     &                0.25_r8*(ubar(i  ,j,krhs)*ubar(i  ,j,krhs)+       &
     &                         ubar(i+1,j,krhs)*ubar(i+1,j,krhs)+       &
     &                         vbar(i,j  ,krhs)*vbar(i,j  ,krhs)+       &
     &                         vbar(i,j+1,krhs)*vbar(i,j+1,krhs))
.....
          avgke=0.5_r8*avgke/volume
          avgpe=avgpe/volume
          avgkp=avgke+avgpe

To be corret, one must either use 0.25_r8*(ubar(i ,j,krhs)*... and avgke=avgke/volume Or use 0.5_r8*(ubar(i ,j,krhs)*... and avgke=0.5_r8*avgke/volume. This problem has been discussed some days ago.[/code]

arango · #5 Unread post by **arango** » Sun Nov 19, 2006 12:37 am

Yes

, I don't recall how the density was added to the kinitic energy. We need to use the 0.25 factor when computing ke2d. We also just need to divide the averaged kinetic enegy by just the volume, avgke=avgke/volume.

This is a diagnostic quantity and it is not used in ROMS. It just printed to standard output. I corrected this bug and added the corrections tar file for version 2.2. Thank you for finding and reporting this bug.