[MITgcm-devel] Inverted echo sounders
Matthew Mazloff
mmazloff at ucsd.edu
Wed Jan 19 16:15:01 EST 2011
Hi Martin,
The important pieces of code to check in are
ctrl_getobcse.F ctrl_getobcsn.F ctrl_getobcss.F ctrl_getobcsw.F
ctrl_init_obcs_variables.F
These all run properly with my setup. Model runs will be unchanged
unless ALLOW_OBCS_CONTROL_MODES is defined.
It would be nice to also add the simple additions in
ecco_cost_final.F too
Thanks!
Matt
On Jan 19, 2011, at 12:37 PM, Matthew Mazloff wrote:
> Hi Martin,
>
> Thanks for checking that in!
>
> So I'm merging your changes and also ensuring backward compatibility
> with ALLOW_CTRL_OBCS_BALANCE (even though I'm not sure the status of
> that piece of code).
>
> I'll check in my updates after they've been checked...
>
> -Matt
>
> ps> you put an #if instead of an #ifdef on line 94 of
> ctrl_get_gen_rec.F
>
>
> On Jan 19, 2011, at 5:59 AM, Martin Losch wrote:
>
>> Hi Matt,
>> as you may have seen I have checked in the cost_obcs*.F files (I
>> modified some comments to be a little more in line with general
>> standards, I hope) from your directory. These were exactly the same
>> files that I have been using for a while now, and never got around
>> to checking in.
>>
>> This moring, I also checked in changes to ctrl_getobcs*.F that
>> should be compatible with what you have in MITgcm_contrib/SOSE/
>> BoxAdj/code_ad/ (basically got rid of exf_swapfields), but you need
>> to merge them.
>>
>> Martin
>>
>> On Jan 18, 2011, at 8:50 PM, Matthew Mazloff wrote:
>>
>>> Hi Patrick,
>>>
>>> This code:
>>> /u/gcmpack/MITgcm_contrib/SOSE/BoxAdj
>>> is compatible with the latest code.
>>>
>>> You can compile it and run it with the files in 'input'. (gendata
>>> is included too)
>>>
>>> The only cost term right now is cost_theta. To actually see the
>>> true benefit of this code one must also enable ssh cost
>>> constraints, which then make the barotropic sensitivity swamp the
>>> baroclinic (cost_theta) sensitivity without this code. But I
>>> wanted to keep the setup simple for now...let me know if you want
>>> me to add this though
>>>
>>> Anyway, it works. The only place the code could use the help of an
>>> expert is for subroutine
>>> ctrl_init_obcs_variables.F
>>> where the reading in of the modes is unrefined...perhaps you can
>>> figure out a cleaner way to do it
>>>
>>> Let me know what you think -- it would be great to have this in the
>>> main branch
>>>
>>> Thanks
>>> Matt
>>>
>>>
>>> On Jan 18, 2011, at 4:47 AM, Patrick Heimbach wrote:
>>>
>>>> Hi Matt,
>>>>
>>>> thanks for the appended code.
>>>> The most promising way to get this checked in cleanly would be to
>>>> 1. initially put the relevant code and code modifs in
>>>> MITgcm_contrib/
>>>> 2. provide an (as simple as possible) working reference setup as a
>>>> test config.
>>>>
>>>> Do you think you could do these two parts?
>>>> It would then be much simpler to merge it into the main branch
>>>> (works often, but not always...).
>>>>
>>>> Cheers
>>>> -Patrick
>>>>
>>>> On Jan 17, 2011, at 12:34 PM, Matthew Mazloff wrote:
>>>>
>>>>> Hi Martin,
>>>>>
>>>>> That, but also the ability to control the obcs using a modal
>>>>> decomposition. This is important when optimizing to both SSH and
>>>>> in situ obs as barotropic sensitivities from SSH cost overwhelms
>>>>> baroclinic sensitivities. Thus the only way to optimize both is
>>>>> to separate the components and pre-condition accordingly.
>>>>>
>>>>> The code changes involve
>>>>> -- adding the mode matrix to ctrl.h
>>>>> -- reading in the mode matrix in ctrl_init_obcs_variables.F
>>>>> -- changing ctrl_getobcs[n,s,w,e].F to do the decomposition
>>>>>
>>>>> Its fairly simple...it would be great if you wanted to check it
>>>>> in.
>>>>>
>>>>> The only issue is one likely would want to fix the reading in of
>>>>> the modes to be more consistent with other MITgcm I/O practices.
>>>>> I will attach this piece of the code below. Also, in this part
>>>>> of the code I put the only "Read Me" aspect...so a lot of
>>>>> (hopefully useful) comments
>>>>>
>>>>> Let me know what you think,
>>>>> Matt
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> This is starting at line ~131 of ctrl_init_obcs_variables.F
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> #ifdef ALLOW_OBCS_CONTROL
>>>>> cih Get matrices for reconstruction from barotropic-barclinic
>>>>> modes
>>>>> CMM( Likely want to define filename=baro_invmodes.bin at runtime
>>>>> c using data.ctrl...hardcoded with no checks for now...)
>>>>> c if baro_invmodes.bin is given then will use modes
>>>>> c otherwise z and k space equal -- no modes
>>>>> c really need more checks on this though!!!!
>>>>>
>>>>> do k = 1,nr
>>>>> do j = 1,nr
>>>>> do i = 1,nr
>>>>> modesv(i,j,k) = 1. _d 0
>>>>> enddo
>>>>> enddo
>>>>> enddo
>>>>> inquire( file='baro_invmodes.bin', exist=exst )
>>>>> IF ( exst ) THEN
>>>>>
>>>>> CMM -- i want to use mdsreadvector but cannot as it assumes some
>>>>> c grid info and will increment record based on processor..
>>>>> c so for now use this -- it works -- and fix in future....
>>>>>
>>>>> open(117, file='baro_invmodes.bin', access='direct',
>>>>> & recl=2*WORDLENGTH*Nr*Nr, status='old')
>>>>> do nz = 1,Nr
>>>>> read(117,rec=nz) ((modesv(k,nk,nz), k=1,Nr), nk=1,Nr)
>>>>> end do
>>>>>
>>>>> CMM AND HERE IS README PART OF CODE:
>>>>> CMM modesv should be orthonormal modes -- ideally solution of
>>>>> c planetary vertical mode equation using model mean dRho/dz
>>>>> c modesv is nr X nr X nr
>>>>> c dim one is z-space
>>>>> c dim two is mode space
>>>>> c dim three is the total depth for which this set of modes
>>>>> applies
>>>>> c so for example modesv(:,2,nr) will be the second mode
>>>>> c in z-space for the full model depth
>>>>> c It is important to note that the modes must be orthogonal
>>>>> c when weighted by dz!
>>>>> c i.e. if f_i(z) = mode i, sum_j(f_i(z_j)*f_j(z_j)*dz_j =
>>>>> delta_ij
>>>>> c first mode should also be constant in depth...barotropic
>>>>> c here is matlab code to ensure modes are orthonormal
>>>>> c ________________________________________________________
>>>>> C for k = 2:NZ;
>>>>> C iz = 1:k;% iz = vector of depth levels
>>>>> C NZ-k % print out a progress number
>>>>> C % have to weight by dz! Make a vector of fractional dz's
>>>>> C zwt = DRF(iz)/sum(DRF(iz),1);
>>>>> C %ensure first mode is barotropic (constant in depth)
>>>>> C avm1=sum(modesv(iz,1,k).*zwt,1);
>>>>> C modesv(iz,1,k)=avm1;
>>>>> C %make all modes orthogonal weighted by delta z
>>>>> C % use gram-schmidt leaving first one the same
>>>>> C for mds = 1:k-1
>>>>> C R = sum((modesv(iz,mds,k).^2).*zwt,1);
>>>>> C R2 = (modesv(iz,mds,k).*zwt)'*modesv(iz,mds+1:k,k)/R;
>>>>> C modesv(iz,mds+1:k,k) = modesv(iz,mds+1:k,k) - ...
>>>>> C modesv(iz,mds,k)*R2;
>>>>> C end
>>>>> C %All now orthognal, now normalize
>>>>> C for mds = 1:k
>>>>> C R = sqrt(sum((modesv(iz,mds,k).^2).*zwt,1));
>>>>> C if R < 1e-8
>>>>> C fprintf('Small R!! for mds = %2i and k = %2i
>>>>> \n',mds,k);
>>>>> C R = inf;
>>>>> C end
>>>>> C modesv(iz,mds,k) = modesv(iz,mds,k)./R;
>>>>> C end
>>>>> C end;% end of loop over depth level k
>>>>> C fid = fopen('baro_invmodes.bin','w','b');
>>>>> C fwrite(fid,modesv,'single');%or double depending on control
>>>>> prec
>>>>> C fclose(fid)
>>>>> C if 1 %plot first 5 modes for deepest case
>>>>> C figure
>>>>> C clf;plot(modesv(:,[1:5],NZ),RC(:));
>>>>> C title('output modes at deepest point')
>>>>> C end
>>>>> C if 1 % test orthogonality
>>>>> C % do whole matrix (need to include dz!)
>>>>> C k = NZ;
>>>>> C cmm = (squeeze(modesv(iz,iz,k)).*repmat(zwt,[1 k]))'* ...
>>>>> C squeeze(modesv(iz,iz,k));
>>>>> C figure;imagesc(cmm);colorbar
>>>>> C title([num2str(k) ' mode orthogonality min, max diag:
>>>>> ' ...
>>>>> C num2str(min(diag(cmm))) ', ' ...
>>>>> C num2str(max(diag(cmm)))])
>>>>> C end
>>>>> C save baro_invmodes modesv
>>>>> c ________________________________________________________
>>>>>
>>>>> #endif
>>>>>
>>>>>
>>>>> return
>>>>> end
>>>>>
>>>>>
>>>>>
>>>>> ~
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On Jan 17, 2011, at 12:29 AM, Martin Losch wrote:
>>>>>
>>>>>> Matt,
>>>>>>
>>>>>> I have obcs code to check in: regularization terms for the ECCO
>>>>>> cost function, which I got from you (ecco_cost_obcs[n,s,w,e].F).
>>>>>> Is that the change that you would like to see in the repository?
>>>>>>
>>>>>> Martin
>>>>>>
>>>>>> On Jan 15, 2011, at 12:27 AM, Matthew Mazloff wrote:
>>>>>>
>>>>>>> ps> would like to get obcs mode code checked in too...
>>>>>>
>>>>>>
>>>>>> _______________________________________________
>>>>>> 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
>>>>
>>>> ---
>>>> Patrick Heimbach | heimbach at mit.edu | http://www.mit.edu/~heimbach
>>>> MIT | EAPS 54-1518 | 77 Massachusetts Ave | Cambridge MA 02139 USA
>>>> FON +1-617-253-5259 | FAX +1-617-253-4464 | SKYPE patrick.heimbach
>>>>
>>>>
>>>>
>>>> _______________________________________________
>>>> 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