[MITgcm-devel] woes with exf interpolating

Matthew Mazloff mmazloff at ucsd.edu
Mon Apr 5 19:22:12 EDT 2010


Hi JMC,

Thanks for your help!  I had some trouble checking in SOSE but will  
try to revisit later this week.

As for your questions, I do compile with package debug.  I do not  
TARGET_NEC_SX , could that help?  And yeah,   xgOrigin = 0.0,

I have not checked regarding optimization -- can I just NOOPTFILE  
exf_interp.F or are there other routines I should include?

And I can try your fix to see if anything changes...

Thanks!
-Matt








On Apr 5, 2010, at 11:23 AM, Jean-Michel Campin wrote:

> Hi Matt,
>
> I did take a look at exf_interp.F, and did not see any obvious
> reason for a problem.
> One thing that would help is to have your set-up checked-in,
> (I remember Patrick creating this MITgcm_contrib/SOSE 2 months ago,
> but it's still almost empty).
> but in the mean time, could be useful to know the status of few
> CPP options (TARGET_NEC_SX, ALLOW_DEBUG ...), and I assume
> you are using the default xgOrigin=0 (in "data") ?
>
> Compiling the pkg "debug" would turn on '#define ALLOW_DEBUG'
> and would check for valid range of long & lat in S/R exf_interp,
> which could be useful.
>
> Also, did you check that with no compiler optimization you get
> the same thing and same problem ? (there are few loop sequence,
> with j loop inside an outer i loop that might be changed by the
> compiler).
>
> There is also something you could try in exf_interp.F,
> to put xG in interval [lon0-long_inc*0.5 , lon0-long_inc*0.5 + 360 [
> instead of in [ lon_0 , lon_0+360 [
> by changing lines 113-115, from:
> 113           xG(i,j,bi,bj) = xG_in(i,j,bi,bj)-lon_0
> 114      &                  + threeSixtyRS*2.
> 115           xG(i,j,bi,bj) = lon_0+mod(xG(i,j,bi,bj),threeSixtyRS)
> to the modified version:
> 113           xG(i,j,bi,bj) = xG_in(i,j,bi,bj)-lon_0+long_inc*0.5 _d 0
> 114      &                  + threeSixtyRS*2.
> 115           xG(i,j,bi,bj) = lon_0-long_inc*0.5 _d 0
> 116      &                  + mod(xG(i,j,bi,bj),threeSixtyRS)
> It will be more "symetric" (although I don't think it should be
> a problem, but I could be wrong).
>
> Cheers,
> Jean-Michel
>
> On Sat, Apr 03, 2010 at 12:22:47PM -0700, Matthew Mazloff wrote:
>> Thanks Dimitris!  I'll look for other explanations
>>
>> appreciate it,
>> -Matt
>>
>>
>> On Apr 3, 2010, at 12:14 PM, Menemenlis, Dimitris (3248) wrote:
>>
>>> Matt, exf_interp assumes periodicity in zonal direction of input  
>>> array
>>> and provides 2 overlaps on either side, which should be sufficient  
>>> for
>>> bicubic, so I don't understand where discontinuity comes from.
>>>
>>> C enlarge boundary
>>>     do j=1,ny_in
>>>      arrayin(0,j)       = arrayin(nx_in,j)
>>>      arrayin(-1,j)      = arrayin(nx_in-1,j)
>>>      arrayin(nx_in+1,j) = arrayin(1,j)
>>>      arrayin(nx_in+2,j) = arrayin(2,j)
>>>     enddo
>>>     do i=-1,nx_in+2
>>>      arrayin(i,0)       = arrayin(i,1)
>>>      arrayin(i,-1)      = arrayin(i,1)
>>>      arrayin(i,ny_in+1) = arrayin(i,ny_in)
>>>      arrayin(i,ny_in+2) = arrayin(i,ny_in)
>>>     enddo
>>>
>>> A simple test would be to add one more overlap on each side of
>>> arrayin, and also to the related x_in and y_in, and see if it fixes
>>> problem.
>>>
>>> Dimitris
>>>
>>> On Apr 3, 2010, at 10:45 AM, Matthew Mazloff wrote:
>>>
>>>> Hi JMC,
>>>>
>>>> I copied it below -- I think its pretty standard.  Is there  
>>>> something
>>>> I'm missing in prescribing the forcing dimensions -- do I need to
>>>> manually provide an overlap, or provide overlap information?
>>>>
>>>> Thanks
>>>> Matt
>>>>
>>>>
>>>>
>>>> login4% more data.exf
>>>> # $Header: /u/gcmpack/MITgcm/verification/lab_sea/input/data.exf,v
>>>> 1.7
>>>> 2004/07/11 06:21:44 dimitri Exp
>>>> $
>>>> #
>>>> # *********************
>>>> # External Forcing Data
>>>> # *********************
>>>> &EXF_NML_01
>>>> #
>>>> useExfYearlyFields = .TRUE.,
>>>> useExfCheckRange = .FALSE.,
>>>> exf_iprec         = 32,
>>>> exf_yftype        = 'RL',
>>>> useRelativeWind   = .TRUE.,
>>>> /
>>>> #
>>>> &EXF_NML_02
>>>> #
>>>> atempstartdate1   = 20050101,
>>>> atempstartdate2   = 00000,
>>>> atempperiod       = 21600.0,
>>>> #
>>>> aqhstartdate1     = 20050101,
>>>> aqhstartdate2     = 00000,
>>>> aqhperiod         = 21600.0,
>>>> #
>>>> precipstartdate1  = 20050101,
>>>> precipstartdate2  = 00000,
>>>> precipperiod      = 21600.0,
>>>> #
>>>> uwindstartdate1   = 20050101,
>>>> uwindstartdate2   = 00000,
>>>> uwindperiod       = 21600.0,
>>>> #
>>>> vwindstartdate1   = 20050101,
>>>> vwindstartdate2   = 00000,
>>>> vwindperiod       = 21600.0,
>>>> #
>>>> runoffstartdate1  = 20050101,
>>>> runoffstartdate2  = 00000,
>>>> runoffperiod      = 0.0,
>>>> #
>>>> swdownstartdate1  = 20050101,
>>>> swdownstartdate2  = 00000,
>>>> swdownperiod      = 21600.0,
>>>> #
>>>> lwdownstartdate1  = 20050101,
>>>> lwdownstartdate2  = 00000,
>>>> lwdownperiod      = 21600.0,
>>>> #
>>>> atempfile         = 'SO6input/NCEP_R6_tmp2m_degC_SO',
>>>> aqhfile           = 'SO6input/NCEP_R6_spfh2m_SO',
>>>> uwindfile         = 'SO6input/NCEP_R6_u10m_SO',
>>>> vwindfile         = 'SO6input/NCEP_R6_v10m_SO',
>>>> precipfile        = 'SO6input/NCEP_R6_rain_SO_DT',
>>>> runoffFile        = 'SO6input/run-off_SO.6',
>>>> lwdownfile        = 'SO6input/NCEP_R6_dlw_SO',
>>>> swdownfile        = 'SO6input/NCEP_R6_dsw_SO',
>>>> #
>>>> /
>>>> #
>>>> &EXF_NML_03
>>>> #
>>>> exf_inscal_uwind = 1.,
>>>> exf_inscal_vwind = 1.,
>>>> exf_inscal_runoff = 3.170979198E-8,
>>>> exf_offset_atemp = 273.15,
>>>> #
>>>> /
>>>> #
>>>> &EXF_NML_04
>>>> #
>>>> atemp_lon0       = 0.5D0,
>>>> atemp_lon_inc     = 1.D0,
>>>> atemp_lat0        = -79.5D0,
>>>> atemp_lat_inc     = 60*1.D0,
>>>> atemp_nlon        = 360,
>>>> atemp_nlat        = 60,
>>>> #
>>>> aqh_lon0          = 0.5D0,
>>>> aqh_lon_inc       = 1.D0,
>>>> aqh_lat0          = -79.5D0,
>>>> aqh_lat_inc       = 60*1.D0,
>>>> aqh_nlon          = 360,
>>>> aqh_nlat          = 60,
>>>> #
>>>> precip_lon0       = 0.5D0,
>>>> precip_lon_inc    = 1.D0,
>>>> precip_lat0       = -79.5D0,
>>>> precip_lat_inc    = 60*1.D0,
>>>> precip_nlon       = 360,
>>>> precip_nlat       = 60,
>>>> #
>>>> uwind_lon0         = 0.5D0,
>>>> uwind_lon_inc      = 1.D0,
>>>> uwind_lat0         = -79.5D0,
>>>> uwind_lat_inc      = 60*1.D0,
>>>> uwind_nlon         = 360,
>>>> uwind_nlat         = 60,
>>>> #
>>>> vwind_lon0         = 0.5D0,
>>>> vwind_lon_inc      = 1.D0,
>>>> vwind_lat0         = -79.5D0,
>>>> vwind_lat_inc      = 60*1.D0,
>>>> vwind_nlon         = 360,
>>>> vwind_nlat         = 60,
>>>> #
>>>> swdown_lon0          = 0.5D0,
>>>> swdown_lon_inc       = 1.D0,
>>>> swdown_lat0          = -79.5D0,
>>>> swdown_lat_inc       = 60*1.D0,
>>>> swdown_nlon          = 360,
>>>> swdown_nlat          = 60,
>>>> #
>>>> lwdown_lon0          = 0.5D0,
>>>> lwdown_lon_inc       = 1.D0,
>>>> lwdown_lat0          = -79.5D0,
>>>> lwdown_lat_inc       = 60*1.D0,
>>>> lwdown_nlon          = 360,
>>>> lwdown_nlat          = 60,
>>>> #
>>>> /
>>>> #
>>>> &EXF_NML_OBCS
>>>> #
>>>> obcsNstartdate1     = 20050101,
>>>> obcsNstartdate2     = 00000,
>>>> obcsNperiod         = 2635200.0,
>>>> #
>>>> /
>>>>
>>>>
>>>>
>>>> On Apr 3, 2010, at 7:07 AM, Jean-Michel Campin wrote:
>>>>
>>>>> Hi Matt,
>>>>>
>>>>> Could you send also the "data.exf" you are using ? or could I find
>>>>> it
>>>>> somewhere ?
>>>>> Thanks,
>>>>> Jean-Michel
>>>>>
>>>>> On Fri, Apr 02, 2010 at 03:17:51PM -0700, Matthew Mazloff wrote:
>>>>>> Hello,
>>>>>>
>>>>>> A problem has been found (thanks Fabien!) regarding a  
>>>>>> discontinuity
>>>>>> in
>>>>>> the wind fields around longitude 0/360.
>>>>>>
>>>>>> My set-up is 1/6 of a degree and I use NCEP on a 1 degree grid  
>>>>>> and
>>>>>> use
>>>>>> exf interpolation.  Within ~1 degree east and west of longitude 0
>>>>>> there
>>>>>> is a discontinuity in the wind field -- it is small but quite
>>>>>> significant in the x-derivative.  My guess is that the
>>>>>> interpolation
>>>>>> routine is not properly accounting for periodicity -- and instead
>>>>>> perhaps incorporating zeros.
>>>>>>
>>>>>> Some important stuff -- I am using NCEP winds -- not stress --
>>>>>> and I
>>>>>> have switched to interp_method=2 instead of the hard-coded
>>>>>> interp_method=1.  This is because with interp_method 1 (linear)
>>>>>> there
>>>>>> are discontinuities in the wind-stress curl everywhere.  With
>>>>>> bicubic
>>>>>> interp the curl is smooth -- with this exception around  
>>>>>> longitude 0
>>>>>>
>>>>>> I was hoping someone with more expertise could shed light on this
>>>>>> before
>>>>>> I go digging in the code
>>>>>>
>>>>>> Thanks!
>>>>>> -Matt
>>>>>>
>>>>>>
>>>>>> _______________________________________________
>>>>>> MITgcm-devel mailing list
>>>>>> MITgcm-devel at mitgcm.org
>>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>>>>>
>>>>> _______________________________________________
>>>>> MITgcm-devel mailing list
>>>>> MITgcm-devel at mitgcm.org
>>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>>>>
>>>>
>>>> _______________________________________________
>>>> MITgcm-devel mailing list
>>>> MITgcm-devel at mitgcm.org
>>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>>>
>>>
>>> _______________________________________________
>>> MITgcm-devel mailing list
>>> MITgcm-devel at mitgcm.org
>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>>
>>
>> _______________________________________________
>> MITgcm-devel mailing list
>> MITgcm-devel at mitgcm.org
>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel




More information about the MITgcm-devel mailing list