Changeset 8936


Ignore:
Timestamp:
07/12/11 17:47:39 (14 years ago)
Author:
Eric.Larour
Message:

Implemented Mat Array for qmu.

Location:
issm/trunk/src
Files:
14 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

    r8934 r8936  
    404404        QmuNumberOfResponsesEnum,
    405405        QmuSaveFemmodelEnum,
     406        PartEnum,
    406407        IndexedEnum,
    407408        RegularEnum,
    … …  
    409410        ScaledEnum,
    410411        NodalEnum,
    411         ResponseDescriptorsEnum,
     412        ResponsedescriptorsEnum,
    412413        PetscRcEnum,
    413414        SparsityEnum,
    414415        TolXEnum, //TO BE DELETED
    415         VariableDescriptorsEnum,
     416        VariabledescriptorsEnum,
    416417        VerboseEnum,
    417418        WaitOnLockEnum, //TO BE DELETED
    … …  
    537538        WaitonlockEnum,
    538539        NumberOfVariablesEnum,
    539         NumvariabledescriptorsEnum,
    540540        NumberOfResponsesEnum,
    541         NumresponsedescriptorsEnum,
    542541        NpartEnum,
    543542        QmuMassFluxNumProfilesEnum

  • issm/trunk/src/c/io/Disk/IoModelFetchData.cpp

    r8926 r8936  
    215215        *pinteger=integer;
    216216
     217}
     218/*}}}*/
     219/*FUNCTION IoModelFetchData(double*** pmatrices,int** pmdims,int** pndims, int* pM,FILE* model_handle,int data_enum){{{1*/
     220void  IoModelFetchData(double*** pmatrices,int** pmdims,int** pndims, int* pnumrecords,FILE* model_handle,int data_enum){
     221
     222        int i;
     223
     224        extern int my_rank;
     225        extern int num_procs;
     226
     227        /*output: */
     228        double** matrices=NULL;
     229        int*     mdims=NULL;
     230        int*     ndims=NULL;
     231        int      numrecords=0;
     232
     233        /*intermediary: */
     234        int M,N;
     235        double* matrix=NULL;
     236       
     237        FILE* fid=NULL;
     238       
     239        /*Set file pointer to beginning of the data: */
     240        fid=SetFilePointerToData(model_handle,data_enum);
     241       
     242        /*Now fetch: */
     243        if(my_rank==0){ 
     244                if(fread(&numrecords,sizeof(int),1,fid)!=1) _error_("could not read number of records in matrix array ");
     245        }
     246        MPI_Bcast(&numrecords,1,MPI_INT,0,MPI_COMM_WORLD);
     247
     248        if(numrecords){
     249
     250                /*Allocate matrices :*/
     251                matrices=(double**)xmalloc(numrecords*sizeof(double*));
     252                mdims=(int*)xmalloc(numrecords*sizeof(int));
     253                ndims=(int*)xmalloc(numrecords*sizeof(int));
     254
     255                for(i=0;i<numrecords;i++){
     256                        matrices[i]=NULL;
     257                        mdims[i]=NULL;
     258                        ndims[i]=NULL;
     259                }
     260
     261                /*Loop through records and fetch matrix: */
     262                for(i=0;i<numrecords;i++){
     263
     264                        if(my_rank==0){ 
     265                                if(fread(&M,sizeof(int),1,fid)!=1) _error_("%s%i%s","could not read number of rows in ",i,"th matrix of matrix array");
     266                        }
     267                        MPI_Bcast(&M,1,MPI_INT,0,MPI_COMM_WORLD);
     268
     269                        if(my_rank==0){ 
     270                                if(fread(&N,sizeof(int),1,fid)!=1) _error_("%s%i%s","could not read number of columns in ",i,"th matrix of matrix array");
     271                        }
     272                        MPI_Bcast(&N,1,MPI_INT,0,MPI_COMM_WORLD);
     273
     274                        /*Now allocate matrix: */
     275                        if(M*N){
     276                                matrix=(double*)xmalloc(M*N*sizeof(double));
     277
     278                                /*Read matrix on node 0, then broadcast: */
     279                                if(my_rank==0){ 
     280                                        if(fread(matrix,M*N*sizeof(double),1,fid)!=1) _error_("could not read matrix ");
     281                                }
     282
     283                                MPI_Bcast(matrix,M*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
     284                        }
     285
     286                        /*Assign: */
     287                        matrices[i]=matrix;
     288                        mdims[i]=M;
     289                        ndims[i]=N;
     290                }
     291        }
     292
     293        /*Assign output pointers: */
     294        *pmatrices=matrices;
     295        *pmdims=mdims;
     296        *pndims=ndims;
     297        *pnumrecords=numrecords;
    217298}
    218299/*}}}*/

  • issm/trunk/src/c/io/Disk/diskio.h

    r8926 r8936  
    2121void  IoModelFetchData(int* pinteger,FILE* model_handle,int data_enum);
    2222void  IoModelFetchData(char*** pstrings,int* pnumstrings,FILE* model_handle,int data_enum);
     23void  IoModelFetchData(double*** pmatrices,int** pmdims,int** pndims, int* pnumrecords,FILE* model_handle,int data_enum);
    2324
    2425FILE* SetFilePointerToData(FILE* model_handle,int data_enum);

  • issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp

    r8934 r8936  
    354354                case QmuNumberOfResponsesEnum : return "QmuNumberOfResponses";
    355355                case QmuSaveFemmodelEnum : return "QmuSaveFemmodel";
     356                case PartEnum : return "Part";
    356357                case IndexedEnum : return "Indexed";
    357358                case RegularEnum : return "Regular";
    … …  
    359360                case ScaledEnum : return "Scaled";
    360361                case NodalEnum : return "Nodal";
    361                 case ResponseDescriptorsEnum : return "ResponseDescriptors";
     362                case ResponsedescriptorsEnum : return "Responsedescriptors";
    362363                case PetscRcEnum : return "PetscRc";
    363364                case SparsityEnum : return "Sparsity";
    364365                case TolXEnum : return "TolX";
    365                 case VariableDescriptorsEnum : return "VariableDescriptors";
     366                case VariabledescriptorsEnum : return "VariableDescriptors";
    366367                case VerboseEnum : return "Verbose";
    367368                case WaitOnLockEnum : return "WaitOnLock";
    … …  
    478479                case WaitonlockEnum : return "Waitonlock";
    479480                case NumberOfVariablesEnum : return "NumberOfVariables";
    480                 case NumvariabledescriptorsEnum : return "Numvariabledescriptors";
    481481                case NumberOfResponsesEnum : return "NumberOfResponses";
    482                 case NumresponsedescriptorsEnum : return "Numresponsedescriptors";
    483482                case NpartEnum : return "Npart";
    484483                case QmuMassFluxNumProfilesEnum : return "QmuMassFluxNumProfiles";

  • issm/trunk/src/c/modules/ModelProcessorx/Qmu/CreateParametersQmu.cpp

    r8926 r8936  
    2525
    2626        char**  responsedescriptors=NULL;
     27        int     numresponsedescriptors;
    2728        char**  variabledescriptors=NULL;
     29        int     numvariabledescriptors;
    2830        char*   descriptor=NULL;
    2931        double* dakota_parameter=NULL;
    … …  
    3941        /*parameters for mass flux: */
    4042        bool     qmu_mass_flux_present=false;
     43        int qmu_mass_flux_num_profiles=0;
    4144        double** array=NULL;
     45        int*     mdims_array=NULL;
     46        int*     ndims_array=NULL;
     47       
    4248        double*  matrix=NULL;
    4349        double*  temp_matrix=NULL;
    44         int*     mdims_array=NULL;
    45         int*     ndims_array=NULL;
    4650        int      M;
    4751        int      temp_m,temp_n;
    … …  
    7074                parameters->AddObject(new   StringParam(QmuErrNameEnum,qmuerrname));
    7175                /*}}}*/
    72                 /*Deal with variable descriptors: {{{1*/
    73                 variabledescriptors=(char**)xmalloc(iomodel->numvariabledescriptors*sizeof(char*));
    74 
    75                 /*Fetch descriptors: logic varies if we are running parallel or serial. In parallel, qmumarshall
    76                  * took care of marshalling all the variable descriptors, so it's easy. In serial mode,
    77                  * the variables are in md.variables(md.ivar), as a strucuture: */
    78 
    79                 for(i=0;i<iomodel->numvariabledescriptors;i++){
    80                         sprintf(tag,"%s%i","variabledescriptor",i+1);
    81                         _error_("not yet");
    82                         //IoModelFetchData(&descriptor,iomodel_handle,tag);
    83                         variabledescriptors[i]=descriptor;
    84                 }
     76                /*Fetch variable descriptors: {{{1*/
     77                IoModelFetchData(&variabledescriptors,&numvariabledescriptors,iomodel_handle,VariabledescriptorsEnum);
    8578
    8679                /*Ok, we have all the variable descriptors. Build a parameter with it: */
    87                 parameters->AddObject(new StringArrayParam(VariableDescriptorsEnum,variabledescriptors,iomodel->numvariabledescriptors));
    88 
    89                 /*}}}*/
    90                 /*Deal with response descriptors: {{{1*/
    91                 responsedescriptors=(char**)xmalloc(iomodel->numresponsedescriptors*sizeof(char*));
    92 
    93                 /*Fetch descriptors: */
    94                 for(i=0;i<iomodel->numresponsedescriptors;i++){
    95                         sprintf(tag,"%s%i","responsedescriptor",i+1);
    96                         _error_("not yet");
    97                         //IoModelFetchData(&descriptor,iomodel_handle,tag);
    98                         responsedescriptors[i]=descriptor;
    99                 }
     80                parameters->AddObject(new StringArrayParam(VariabledescriptorsEnum,variabledescriptors,numvariabledescriptors));
     81
     82                /*}}}*/
     83                /*Fetch response descriptors: {{{1*/
     84                IoModelFetchData(&responsedescriptors,&numresponsedescriptors,iomodel_handle,ResponsedescriptorsEnum);
    10085
    10186                /*Ok, we have all the response descriptors. Build a parameter with it: */
    102                 parameters->AddObject(new StringArrayParam(ResponseDescriptorsEnum,responsedescriptors,iomodel->numresponsedescriptors));
     87                parameters->AddObject(new StringArrayParam(ResponsedescriptorsEnum,responsedescriptors,numresponsedescriptors));
    10388                parameters->AddObject(new    IntParam(QmuNumberOfResponsesEnum,iomodel->numberofresponses));
    10489                /*}}}*/
    … …  
    10792                parameters->AddObject(new    IntParam(QmuNPartEnum,iomodel->qmu_npart));
    10893               
    109                 //IoModelFetchData(&dpart,NULL,NULL,iomodel_handle,PartEnum);
    110                 _error_("not yet");
     94                IoModelFetchData(&dpart,NULL,NULL,iomodel_handle,PartEnum);
    11195
    11296                if(!dpart){
    … …  
    122106                /*Deal with data needed because of qmu variables: {{{1*/
    123107               
    124                 for(i=0;i<iomodel->numvariabledescriptors;i++){
     108                for(i=0;i<numvariabledescriptors;i++){
    125109
    126110                        if (strncmp(variabledescriptors[i],"scaled_",7)==0){
    … …  
    129113                               
    130114                                /*Recover data: */
    131                                 //IoModelFetchData(&dakota_parameter,NULL,NULL,iomodel_handle,StringToEnumx(tag));
    132                                 _error_("not yet");
     115                                IoModelFetchData(&dakota_parameter,NULL,NULL,iomodel_handle,StringToEnumx(tag));
    133116
    134117                                /*Convert units: */
    … …  
    144127                /*}}}*/
    145128                /*Deal with data needed to compute qmu responses: {{{1*/
    146                 for(i=0;i<iomodel->numresponsedescriptors;i++){
     129                for(i=0;i<numresponsedescriptors;i++){
    147130                       
    148131                        if(strncmp(responsedescriptors[i],"indexed_MassFlux",16)==0){
    … …  
    154137                if(qmu_mass_flux_present){
    155138                       
    156                         /*We need the segments to be able to compute the mass flux. We have as many groups of segments as we have MassFlux
    157                          *responses. Let's build a DoubleMatArrayParam object with the array of segments: */
    158                        
    159                         if(iomodel->qmu_mass_flux_num_profiles==0)_error_(" qmu_mass_flux_num_profiles is 0, when MassFlux computations were requested!");
    160 
    161                         array=(double**)xmalloc(iomodel->qmu_mass_flux_num_profiles*sizeof(double*));
    162                         mdims_array=(int*)xmalloc(iomodel->qmu_mass_flux_num_profiles*sizeof(int));
    163                         ndims_array=(int*)xmalloc(iomodel->qmu_mass_flux_num_profiles*sizeof(int));
    164                         M=iomodel->qmu_mass_flux_num_profiles;
    165 
    166                         for(i=0;i<M;i++){
    167                                 sprintf(tag,"%s%i","qmu_mass_flux_segments",i+1); //names were created using Matlab indexing convention
    168                                 //IoModelFetchData(&temp_matrix,&temp_m,&temp_n,iomodel_handle,tag);
    169                                 _error_("not yet");
    170 
    171                                 /*This temp_matrix represents all the segments, for all elements. On this cpu, we only have a subset of these
    172                                  * elements. Extract matrix, out of temp_matrix, corresponding to only our elements: */
     139                        /*Fetch the mass flux segments necessary to compute the mass fluxes.  Build a DoubleMatArrayParam object out of them: */
     140                        IoModelFetchData(&array,&mdims_array,&ndims_array,&qmu_mass_flux_num_profiles,iomodel_handle,QmuMassFluxSegmentsEnum);
     141                        if(qmu_mass_flux_num_profiles==0)_error_(" qmu_mass_flux_num_profiles is 0, when MassFlux computations were requested!");
     142
     143                        /*Go through segments, and extract those that belong to this cpu: */
     144                        for(i=0;i<qmu_mass_flux_num_profiles;i++){
     145                                temp_matrix=array[i];
     146                                temp_m=mdims_array[i];
     147                                temp_n=ndims_array[i];
     148
    173149                                m=0;
    174150                                for(j=0;j<temp_m;j++){
    … …  
    194170                                ndims_array[i]=5;
    195171
    196                                 /*Free data: */
     172                                /*Free temporary matrix: */
    197173                                xfree((void**)&temp_matrix);
    198174                        }
    199175
    200176                        /*Ok, we have an array of segments, different on every cpu. Create a DoubleMatArrayParam object with it: */
    201                         parameters->AddObject(new DoubleMatArrayParam(QmuMassFluxSegmentsEnum,array,M,mdims_array,ndims_array));
     177                        parameters->AddObject(new DoubleMatArrayParam(QmuMassFluxSegmentsEnum,array,qmu_mass_flux_num_profiles,mdims_array,ndims_array));
    202178
    203179                        /*Free data: */
    … …  
    211187                }
    212188                /*}}}*/
    213 
    214189                parameters->AddObject(new   BoolParam(QmuSaveFemmodelEnum,iomodel->qmu_save_femmodel));
    215 
    216190                /*Free data: {{{1*/
    217                 for(i=0;i<iomodel->numresponsedescriptors;i++){
     191                for(i=0;i<numresponsedescriptors;i++){
    218192                        descriptor=responsedescriptors[i];
    219193                        xfree((void**)&descriptor);
    … …  
    221195                xfree((void**)&responsedescriptors);
    222196
    223                 for(i=0;i<iomodel->numvariabledescriptors;i++){
     197                for(i=0;i<numvariabledescriptors;i++){
    224198                        descriptor=variabledescriptors[i];
    225199                        xfree((void**)&descriptor);

  • issm/trunk/src/c/modules/Qmux/SpawnCoreParallel.cpp

    r8224 r8936  
    5555       
    5656        /*retrieve parameters: */
    57         femmodel->parameters->FindParam(&responses_descriptors,&numresponsedescriptors,ResponseDescriptorsEnum);
     57        femmodel->parameters->FindParam(&responses_descriptors,&numresponsedescriptors,ResponsedescriptorsEnum);
    5858        femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum);
    5959        femmodel->parameters->FindParam(&control_analysis,ControlAnalysisEnum);

  • issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp

    r8934 r8936  
    352352        else if (strcmp(name,"QmuNumberOfResponses")==0) return QmuNumberOfResponsesEnum;
    353353        else if (strcmp(name,"QmuSaveFemmodel")==0) return QmuSaveFemmodelEnum;
     354        else if (strcmp(name,"Part")==0) return PartEnum;
    354355        else if (strcmp(name,"Indexed")==0) return IndexedEnum;
    355356        else if (strcmp(name,"Regular")==0) return RegularEnum;
    … …  
    357358        else if (strcmp(name,"Scaled")==0) return ScaledEnum;
    358359        else if (strcmp(name,"Nodal")==0) return NodalEnum;
    359         else if (strcmp(name,"ResponseDescriptors")==0) return ResponseDescriptorsEnum;
     360        else if (strcmp(name,"Responsedescriptors")==0) return ResponsedescriptorsEnum;
    360361        else if (strcmp(name,"PetscRc")==0) return PetscRcEnum;
    361362        else if (strcmp(name,"Sparsity")==0) return SparsityEnum;
    362363        else if (strcmp(name,"TolX")==0) return TolXEnum;
    363         else if (strcmp(name,"VariableDescriptors")==0) return VariableDescriptorsEnum;
     364        else if (strcmp(name,"Variabledescriptors")==0) return VariabledescriptorsEnum;
    364365        else if (strcmp(name,"Verbose")==0) return VerboseEnum;
    365366        else if (strcmp(name,"WaitOnLock")==0) return WaitOnLockEnum;
    … …  
    476477        else if (strcmp(name,"Waitonlock")==0) return WaitonlockEnum;
    477478        else if (strcmp(name,"NumberOfVariables")==0) return NumberOfVariablesEnum;
    478         else if (strcmp(name,"Numvariabledescriptors")==0) return NumvariabledescriptorsEnum;
    479479        else if (strcmp(name,"NumberOfResponses")==0) return NumberOfResponsesEnum;
    480         else if (strcmp(name,"Numresponsedescriptors")==0) return NumresponsedescriptorsEnum;
    481480        else if (strcmp(name,"Npart")==0) return NpartEnum;
    482481        else if (strcmp(name,"QmuMassFluxNumProfiles")==0) return QmuMassFluxNumProfilesEnum;

  • issm/trunk/src/c/objects/IoModel.cpp

    r8926 r8936  
    233233        if(this->qmu_analysis){
    234234                IoModelFetchData(&this->numberofvariables,iomodel_handle,NumberOfVariablesEnum);
    235                 IoModelFetchData(&this->numvariabledescriptors,iomodel_handle,NumvariabledescriptorsEnum);
    236235                IoModelFetchData(&this->numberofresponses,iomodel_handle,NumberOfResponsesEnum);
    237                 IoModelFetchData(&this->numresponsedescriptors,iomodel_handle,NumresponsedescriptorsEnum);
    238236                IoModelFetchData(&this->qmu_npart,iomodel_handle,NpartEnum);
    239                 IoModelFetchData(&this->qmu_mass_flux_num_profiles,iomodel_handle,QmuMassFluxNumProfilesEnum);
    240237                IoModelFetchData(&this->qmu_save_femmodel,iomodel_handle,QmuSaveFemmodelEnum);
    241238        }

  • issm/trunk/src/c/objects/IoModel.h

    r8926 r8936  
    200200
    201201                /*qmu: */
    202                 int      qmu_mass_flux_num_profiles;
    203202                int      qmu_save_femmodel;
    204203

  • issm/trunk/src/m/classes/model.m

    r8932 r8936  
    332332                 dakotadat={'',false};
    333333                 qmu_analysis={0,true,'Integer'};
    334                  part={[],false};
     334                 part={[],true,'Mat'};
    335335                 npart={0,true,'Integer'};
    336336                 numberofvariables={0,true,'Integer'};
    337337                 numberofresponses={0,true,'Integer'};
    338                  variabledescriptors={NaN,false};
    339                  responsedescriptors={NaN,false};
    340                  numvariabledescriptors={0,true,'Integer'};
    341                  numresponsedescriptors={0,true,'Integer'};
    342                  qmu_mass_flux_num_profiles={0,true,'Integer'};
     338                 variabledescriptors={{},true,'StringArray'};
     339                 responsedescriptors={{},true,'StringArray'};
    343340                 qmu_mass_flux_profile_directory={NaN,false};
    344341                 qmu_mass_flux_profiles={NaN,false};
    345                  qmu_mass_flux_segments={NaN,false};
     342                 qmu_mass_flux_segments={{},true,'MatArray'};
    346343                 qmu_relax={0,false};
    347344                 qmu_save_femmodel={0,true,'Integer'};

  • issm/trunk/src/m/model/WriteData.m

    r8926 r8936  
    4040        %first write length of record
    4141        fwrite(fid,4,'int');  %1 integer
     42
    4243        %now write integer
    4344        fwrite(fid,data,'int');
    … …  
    5051
    5152        %first get length of string array:
    52         num=length(data);
    53 
     53        num=numel(data);
    5454        %now get length of record:
    5555        recordlength=4; %for length of array
    … …  
    7171                fwrite(fid,string,'char');
    7272        end
     73elseif strcmpi(data_type,'MatArray'),
     74
     75        numrecords=numel(data);
     76       
     77        %first get length of record
     78        recordlength=4; %number of records=
     79        for i=1:numrecords,
     80                matrix=data{i};
     81                s=size(matrix);
     82                recordlength=recordlength+4*2+... %row and col of matrix
     83                                         s(1)*s(2)*8; %matrix of doubles
     84        end
     85
     86        %write length of record
     87        fwrite(fid,recordlength,'int');
     88
     89        %write data, first number of records
     90        fwrite(fid,numrecords,'int');
     91       
     92        %write each matrix:
     93        for i=1:numrecords,
     94                matrix=data{i};
     95                s=size(matrix);
     96                fwrite(fid,s(1),'int');
     97                fwrite(fid,s(2),'int');
     98                fwrite(fid,matrix,'double');
     99        end
    73100else
    74101        error('WriteData error message: data type not supported yet!');

  • issm/trunk/src/m/model/solve.m

    r8577 r8936  
    5959marshall(md);
    6060
    61 %add qmu fields to binary file if running qmu analysis
    62 if md.qmu_analysis,
    63         qmumarshall(md,md.variables(options.ivar),md.responses(options.iresp));
    64 end
    6561
    6662%write a template file for issm to use, in parallel

  • issm/trunk/src/m/qmu/preqmu.m

    r6860 r8936  
    103103md.variabledescriptors=variabledescriptors;
    104104md.responsedescriptors=responsedescriptors;
    105 md.numvariabledescriptors=numel(md.variabledescriptors);
    106 md.numresponsedescriptors=numel(md.responsedescriptors);
    107105
    108106%now, we have to provide all the info necessary for the solutions to compute the responses. For ex, if mass_flux

  • issm/trunk/src/m/qmu/process_qmu_response_data.m

    r5487 r8936  
    4343        end
    4444
    45         %ok, process the domains named in qmu_mass_flux_profiles,  to build a list of segments
     45        %ok, process the domains named in qmu_mass_flux_profiles,  to build a list of segments (MatArray)
    4646        md.qmu_mass_flux_segments=cell(num_mass_flux,1);
    4747
    … …  
    5050        end
    5151
    52         %finally, set qmu_mass_flux_num_profiles, for future marshalling
    53         md.qmu_mass_flux_num_profiles=num_mass_flux;
    5452end
Note: See TracChangeset for help on using the changeset viewer.