Potential Energy and Kinteic Energy use different unit?

Bug reports, work arounds and fixes

Moderators: arango, robertson

Post Reply
Message
Author
bzhang
Posts: 11
Joined: Wed Mar 26, 2003 9:25 pm
Location: CICS/ESSIC University of Maryland

Potential Energy and Kinteic Energy use different unit?

#1 Post by bzhang » Thu Sep 28, 2006 6:59 pm

In diag.F

Code: Select all

          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.

User avatar
shchepet
Posts: 185
Joined: Fri Nov 14, 2003 4:57 pm

#2 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.

User avatar
shchepet
Posts: 185
Joined: Fri Nov 14, 2003 4:57 pm

#3 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
Posts: 11
Joined: Wed Mar 26, 2003 9:25 pm
Location: CICS/ESSIC University of Maryland

#4 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]

User avatar
arango
Site Admin
Posts: 1117
Joined: Wed Feb 26, 2003 4:41 pm
Location: IMCS, Rutgers University
Contact:

#5 Post by arango » Sun Nov 19, 2006 12:37 am

Yes :oops:, 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.

Post Reply