Changeset 967


Ignore:
Timestamp:
06/12/09 16:25:23 (16 years ago)
Author:
Mathieu Morlighem
Message:

use process results for any field

Location:
issm/trunk/src/m/solutions/cielo
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/m/solutions/cielo/diagnostic.m

    r963 r967  
    3838       
    3939                %process results
    40                 md.results.diagnostic=processresults(models,results);
     40                if ~isstruct(md.results), md.results=struct(); end
     41                md.results.diagnostic=processresults(results,models,'diagnostic');
    4142        else
    4243                %launch dakota driver for diagnostic core solution

  • issm/trunk/src/m/solutions/cielo/diagnostic_core.m

    r922 r967  
    111111
    112112%load onto results
     113results.step=1;
     114results.time=0;
    113115results.u_g=u_g;
    114116results.p_g=p_g;

  • issm/trunk/src/m/solutions/cielo/processresults.m

    r959 r967  
    1 function newresults=processresults(models,results);
     1function newresults=processresults(results,models,analysis_type);
    22%PROCESSRESULTS - process results from core solution
    33%
    … …  
    66%
    77%   Usage:
    8 %      results=processresults(md,results,m_dh,m_ds,m_dhu)
     8%      results=processresults(results,models,analysis_type)
    99
    10 %recover models
    11 m_dh=models.dh;
    12 m_ds=models.ds;
    13 m_dhu=models.dhu;
     10%initialize output
     11newresults=struct();
    1412
    15 %recover parameters
    16 debug=m_dhu.parameters.debug;
    17 dim=m_dhu.parameters.dim;
    18 ishutter=m_dhu.parameters.ishutter;
    19 ismacayealpattyn=m_dh.parameters.ismacayealpattyn;
    20 isstokes=m_ds.parameters.isstokes;
    21 yts=m_ds.parameters.yts;
     13%recover models first
     14if strcmpi(analysis_type,'diagnostic'),
     15        m_dh=models.dh;
     16        m_ds=models.ds;
     17        m_dhu=models.dhu;
    2218
    23 %recover results
    24 u_g=results.u_g;
    25 p_g=results.p_g;
     19        %some flags needed
     20        dim=m_dhu.parameters.dim;
     21        ishutter=m_dhu.parameters.ishutter;
     22        ismacayealpattyn=m_dh.parameters.ismacayealpattyn;
     23        isstokes=m_ds.parameters.isstokes;
    2624
    27 newresults=struct();
    28 if dim==2,
    29         newresults.step=1;
    30         newresults.time=1;
    31         newresults.vx=u_g(1:2:end)*yts;
    32         newresults.vy=u_g(2:2:end)*yts;
    33         newresults.vel=sqrt(newresults.vx.^2+newresults.vy.^2);
    34         newresults.pressure=p_g;
     25elseif strcmpi(analysis_type,'thermal'),
     26        m_m=models.m;
     27end
    3528
    36 else
    37        
    38         newresults.step=1;
    39         newresults.time=1;
    40         if isstokes,
    41                 newresults.vx=u_g(1:4:end)*yts;
    42                 newresults.vy=u_g(2:4:end)*yts;
    43                 newresults.vz=u_g(3:4:end)*yts;
    44                 newresults.pressure=p_g;
    45         else
    46                 newresults.vx=u_g(1:3:end)*yts;
    47                 newresults.vy=u_g(2:3:end)*yts;
    48                 newresults.vz=u_g(3:3:end)*yts;
    49                 newresults.pressure=p_g;
     29
     30%go through results
     31for i=1:length(results),
     32
     33        %get field names of results
     34        resultsname=fieldnames(results(i));
     35
     36        %go through all field
     37        for j=1:length(resultsname),
     38
     39                if strcmpi(resultsname{j},'u_g'),
     40
     41                        u_g=results(i).u_g;
     42                        yts=m_ds.parameters.yts;
     43
     44                        if dim==2,
     45                                newresults(i).step=results(i).step;
     46                                newresults(i).time=results(i).time;
     47                                newresults(i).vx=u_g(1:2:end)*yts;
     48                                newresults(i).vy=u_g(2:2:end)*yts;
     49                                newresults(i).vel=sqrt(newresults(i).vx.^2+newresults(i).vy.^2);
     50
     51                        else
     52                                newresults(i).step=results(i).step;
     53                                newresults(i).time=results(i).time;
     54                                if isstokes,
     55                                        newresults(i).vx=u_g(1:4:end)*yts;
     56                                        newresults(i).vy=u_g(2:4:end)*yts;
     57                                        newresults(i).vz=u_g(3:4:end)*yts;
     58                                else
     59                                        newresults(i).vx=u_g(1:3:end)*yts;
     60                                        newresults(i).vy=u_g(2:3:end)*yts;
     61                                        newresults(i).vz=u_g(3:3:end)*yts;
     62                                end
     63                                newresults(i).vel=sqrt(newresults(i).vx.^2+newresults(i).vy.^2+newresults(i).vz.^2);
     64                        end
     65
     66                elseif strcmpi(resultsname{j},'p_g'),
     67
     68                        p_g=results(i).p_g;
     69                        newresults(i).step=results(i).step;
     70                        newresults(i).time=results(i).time;
     71                        newresults(i).pressure=p_g;
     72
     73                elseif strcmpi(resultsname{j},'t_g'),
     74
     75                        t_g=results(i).t_g;
     76                        newresults(i).step=results(i).step;
     77                        newresults(i).time=results(i).time;
     78                        newresults(i).temperature=t_g;
     79
     80                elseif strcmpi(resultsname{j},'m_g'),
     81
     82                        m_g=results(i).m_g;
     83                        yts=m_m.parameters.yts;
     84                        newresults(i).step=results(i).step;
     85                        newresults(i).time=results(i).time;
     86                        newresults(i).melting=m_g*yts;
     87
     88                elseif strcmpi(resultsname{j},'h_g'),
     89
     90                        h_g=results(i).h_g;
     91                        newresults(i).step=results(i).step;
     92                        newresults(i).time=results(i).time;
     93                        newresults(i).thickness=h_g;
     94
     95                elseif strcmpi(resultsname{j},'s_g'),
     96
     97                        s_g=results(i).s_g;
     98                        newresults(i).step=results(i).step;
     99                        newresults(i).time=results(i).time;
     100                        newresults(i).surface=s_g;
     101
     102                elseif strcmpi(resultsname{j},'b_g'),
     103
     104                        b_g=results(i).b_g;
     105                        newresults(i).step=results(i).step;
     106                        newresults(i).time=results(i).time;
     107                        newresults(i).bed=b_g;
     108
     109                elseif strcmpi(resultsname{j},'param_g'),
     110
     111                        param_g=results(i).param_g;
     112                        newresults(i).step=results(i).step;
     113                        newresults(i).time=results(i).time;
     114                        newresults(i).parameter=param_g(1:2:end);
     115
     116                else
     117
     118                        newresults(i).step=results(i).step;
     119                        newresults(i).time=results(i).time;
     120                        newresults(i).(resultsname{j})=results(i).(resultsname{j});
     121
     122                end
    50123        end
    51         newresults.vel=sqrt(newresults.vx.^2+newresults.vy.^2+newresults.vz.^2);
    52124end
Note: See TracChangeset for help on using the changeset viewer.