[MITgcm-devel] smoothing package for dimensional controls

Gael Forget gforget at MIT.EDU
Thu May 19 14:05:50 EDT 2011 

Hi Matt,

thanks for accommodating my concern by untying the added option from SMOOTH_CORREL.
Sorry to be giving you such a hard time with this, but I am just trying to avoid additional confusion.

For further clarity, would it be ok to move this stuff to ctrl_set_unpack_xy.F etc. where 
the reverse ALLOW_NONDIMENSIONAL_CONTROL_IO stuff is done. In fact,
I would suggest ALLOW_DIMENSIONAL_CONTROL_IO as the CPP option name
to make it clear that this is the reverse of ALLOW_NONDIMENSIONAL_CONTROL_IO.

So, for example, the added code for ctrl_set_unpack_xy.F would read as
#ifdef ALLOW_DIMENSIONAL_CONTROL_IO
                     if ( lxxadxx ) then
                        globfld3d(i,bi,ip,j,bj,jp,k) =
     &                       globfld3d(i,bi,ip,j,bj,jp,k)
     &                       * sqrt(globfld2d(i,bi,ip,j,bj,jp))
     &                       / forcingPrecond
                     else
                        globfld3d(i,bi,ip,j,bj,jp,k) =
     &                       globfld3d(i,bi,ip,j,bj,jp,k)
     &                       / sqrt(globfld2d(i,bi,ip,j,bj,jp))
     &                       * forcingPrecond
                     endif
#endif
and would be added this at line 217 just before #ifndef ALLOW_SMOOTH_CORREL2D 
so that you can use it regardless of #ifndef/ifdef ALLOW_SMOOTH_CORREL2D.

Makes sense?

Cheers,
Gael

ps:	the reason why I did the smooth_correl2d stuff during ctrl_map rather than ctrl_unpack
	is that the smoother part needs to be seen by taf, and the dimensionalization 
	comes after the smoother (i.e. on the l.h.s.) in the Weaver and Courtier 01. 


On May 18, 2011, at 10:39 PM, Matthew Mazloff wrote:

> Hi Gael,
> 
> OK, I simply took it outside of ALLOW_SMOOTH_CORREL until I can think of a better way to do it
> 
> I switched the CPP flag to NONDIMENSIONALIZE_CONTROL and have this flag only nondim the controls (for later use in smoothing package).
> 
> An issue is the I left smooth_correl2Dw.F the way it was and use it to non-dim 2d ctrls.  Ideally smooth_correl2Dw would be reverted and a separate function from package ctrl would be used to non-dimensionalize...
> 
> Is what I have done ok for now?
> 
> Thanks
> Matt
> 
> 
> 
> 
> 
> 
> 
> 
> On May 18, 2011, at 6:29 PM, Gael Forget wrote:
> 
>> Hi Matt,
>> 
>>> I don't agree that the important smoothing correlation operator of Weaver and Courtier 2001 requires non-dim controls.  Just let S be a diagonal identity matrix instead of the diagonal uncertainty and all is ok...the important part, that the correlation matrix B is represented by a diagonal matrix times a smoother (B^(0.5) = A L^(0.5) W^(-0.5) ) remains...or am I missing something...
>> I agree with your thought exp. The l.h.s. factor S is the one that sets the 
>> error variance of the dimensional dx. You can set it to one, sure.
>> 
>> But I was questioning what you did on the r.h.s. as part of B^(0.5), rather than 
>> what is done on its l.h.s. of it. Now I am unsure (I would need to redo the 
>> maths) but my recollection was that one cannot mess with the r.h.s. if one 
>> wants to ensure that the error variance of dx be the one spec. by S^2.
>> May be I am mistaken though.
>> 
>> Gael
>> 
>>> Regardless, I can move the ALLOW_DIMENSIONAL_CONTROL stuff outside of ALLOW_SMOOTH_CORREL...but would it be ok to leave the smooth_correl2Dw.F the way it is....or do you object to that?  Or am I missing a better way to do this...
>>> 
>>> Sorry for causing trouble -- and perhaps missing something key -- I can call you if you want to discuss the best way to implement this 
>>> 
>>> Thanks
>>> Matt
>>> 
>>> 
>>> On May 18, 2011, at 5:17 PM, Gael Forget wrote:
>>> 
>>>> 
>>>> Hi Matt,
>>>> 
>>>> I guess I was not fast enough to object, but I am not excited by what you did.
>>>> 
>>>> To smooth dimensional controls I would much prefer you added calls to 
>>>> smooth_diff2d / 3d.F outside of the ifdef ALLOW_SMOOTH_CORREL3D/2D.
>>>> 
>>>> To make my reasons more clear I shall summarize the difference between the 
>>>> smoothing and correlation components of pkg/smooth. I reckon I should have 
>>>> documented this stuff in MITgcm-manual already, so I tried to write the 
>>>> following as a draft for the MITgcm-manual.
>>>> 
>>>> --- 
>>>> 
>>>> The very smoothing boils down to a call to smooth_diff3d / 2d.F
>>>> which can of course happen outside of ALLOW_SMOOTH_CORREL3D/2D.
>>>> The smoothing is done by time integrating diffusion equation with diffusion coeffs
>>>> matching a smoothing spatial scale (typically specified through data.smooth).
>>>> 
>>>> The covariance construct of Weaver and Courtier 2001 (se their Eq.52&15)
>>>> is more elaborate than that. In matrix form, it reads as
>>>>    dx = S B^(0.5) du
>>>> where dx is the dimensional control vector adjustment, 
>>>> du is the nondimensional one, S is the diagonal uncertainty 
>>>> matrix (i.e. the scaling factor), and B is the correlation matrix
>>>>    B^(0.5) = A L^(0.5) W^(-0.5)
>>>> where W is the diagonal matrix of grid cell volumes, L is the smoother, 
>>>> and A is the diagonal normalization matrix. A must be computed 
>>>> so that B=B^(T.5) B^(0.5) has ones on the diagonal. Hence
>>>> B does qualify as a correlation matrix, and B^(T.5) S^(2) B^(0.5) does 
>>>> qualify as a covariance matrix. This is what is done within ifdef 
>>>> ALLOW_SMOOTH_CORREL3D / 2D (assuming I did not put bugs or typos).
>>>> 
>>>> ---
>>>> 
>>>> So, coming back to your checkin, Matt,  I would rather not have CPP brackets 
>>>> within ALLOW_SMOOTH_CORREL3D / 2D that would break this defining property of B. 
>>>> Such as adding S^(-1) on the r.h.s. of B^(0.5) …
>>>> 
>>>> Outside of ALLOW_SMOOTH_CORREL3D / 2D such self-consistency does not 
>>>> matter much, and there is a number of preconditioner options (combinations of 
>>>> ALLOW_NONDIMENSIONAL_CONTROL_IO, CTRL_PACK_PRECISE,
>>>> CTRL_DELZNORM, CTRL_SMOOTH, etc?). 
>>>> 
>>>> My guess is that adding some smooth_diffXd.F to that mix would fit your usage. Would it? 
>>>> 
>>>> Cheers,
>>>> Gael
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> <ATT00001.txt>
>>> 
>>> _______________________________________________
>>> MITgcm-devel mailing list
>>> MITgcm-devel at mitgcm.org
>>> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>> 
>> <ATT00001.txt>
> 
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel

More information about the MITgcm-devel mailing list