InputUpdateFromDakotax.h File Reference

header file for updating datasets from inputs More...

#include "../../classes/classes.h"

Go to the source code of this file.

Functions
void	InputUpdateFromDakotax (FemModel *femmodel, double *variables, char **variables_descriptors, int numvariables)

Detailed Description

header file for updating datasets from inputs

Definition in file InputUpdateFromDakotax.h.

Function Documentation

InputUpdateFromDakotax()

void InputUpdateFromDakotax	(	FemModel *	femmodel,
		double *	variables,
		char **	variables_descriptors,
		int	numvariables
	)

Definition at line 12 of file InputUpdateFromDakotax.cpp.

                                                                                                                     {
 
   int     i,j,k,l;
 
   int     numberofvertices;
   int     numberofelements;
   int     nrows;
   int     ncols;
   IssmDouble **variable_partitions         = NULL;
   IssmDouble * variable_partition         = NULL;
   int * variable_partitions_npart         = NULL;
   int * variable_partitions_nt         = NULL;
   int          variable_partitions_num;
   int          npart;
   int          nt;
   int variablecount=0;
 
   double *distributed_values = NULL;
   double *parameter          = NULL;
   char   *descriptor         = NULL;
   char    root[50]; //root name of variable, ex: DragCoefficent, RhoIce, etc ...
 
   if (VerboseQmu())_printf0_("dakota variables updates\n");
 
   /*retrieve parameters: */
   femmodel->parameters->FindParam(&variable_partitions,&variable_partitions_num,NULL,NULL,QmuVariablePartitionsEnum); 
   femmodel->parameters->FindParam(&variable_partitions_npart,NULL,NULL,QmuVariablePartitionsNpartEnum); 
   femmodel->parameters->FindParam(&variable_partitions_nt,NULL,NULL,QmuVariablePartitionsNtEnum); 
   
   numberofvertices=femmodel->vertices->NumberOfVertices();
   numberofelements=femmodel->elements->NumberOfElements();
 
   /*Go through all dakota descriptors, ex: "rho_ice","thermal_conductivity","thickness1","thickness2", etc ..., and 
    * for each descriptor, take the variable value and plug it into the inputs: */
   for(i=0;i<numdakotavariables;i++){ //these are the dakota variables, for all partitions. 
 
      descriptor=variables_descriptors[i];
   
      if (VerboseQmu())_printf0_("   updating variable " << descriptor << "\n");
 
      /*From descriptor, figure out if the variable is scaled, indexed, nodal, or just a simple variable: */
      if (strncmp(descriptor,"scaled_",7)==0){
      
         /*recover partition vector: */
         variable_partition=variable_partitions[variablecount];
         npart=variable_partitions_npart[variablecount];
         nt=variable_partitions_nt[variablecount];
 
         /*Variable is scaled. Determine root name of variable (ex: scaled_DragCoefficient_1 -> DragCoefficient). Allocate distributed_values and fill the 
          * distributed_values with the next npart variables: */
 
         //strcpy(root,strstr(descriptor,"_")+1); *strstr(root,"_")='\0';
         memcpy(root,strstr(descriptor,"_")+1,(strlen(strstr(descriptor,"_")+1)+1)*sizeof(char));
         *strstr(root,"_")='\0';
 
         distributed_values=xNew<double>(npart*nt);
         for(j=0;j<npart*nt;j++){
            distributed_values[j]=variables[i+j];
         }
 
         
         /*Now, pick up the parameter corresponding to root: */
         femmodel->parameters->FindParamInDataset(&parameter,&nrows,&ncols,QmuVariableDescriptorsEnum,StringToEnumx(root));
 
         /*We've got the parameter, we need to update it using the partition vector, and the distributed_values. 
          In addition, the parameter can be either a static or transient (nrows+1) vector. Finally, the partition vectors can include
          -1 (meaning, don't update). */
         
         //_printf_("nrows: " << nrows << " numberofvertices: " << numberofvertices << " numberofelements: " << numberofelements << "\n");
 
         if(ncols!=nt){
            /*we are trying to update a col sized transient input by scaling with a matrix of col size nt. This can only work if nt==1, otherwise, error out: */
            if (nt!=1) _error_("InputUpdateFromDakotax error message: transient input being updated should be the same col size as the number of time step in the qmu variable specificationi");
         }
 
         if(nt==1){
            /*scale all the columns by the same vector:*/
            if (nrows==numberofvertices || nrows==(numberofvertices+1)){
               for(k=0;k<numberofvertices;k++){
                  if (variable_partition[k]==-1)continue;
                  else{
                     for(l=0;l<ncols;l++){
                        *(parameter+ncols*k+l)=*(parameter+ncols*k+l)*distributed_values[(int)variable_partition[k]];
                     }
                  }
               }
            }
            else if (nrows==numberofelements || nrows==(numberofelements+1)){
               for(k=0;k<numberofelements;k++){
                  if (variable_partition[k]==-1)continue;
                  else{
                     for(l=0;l<ncols;l++){
                        *(parameter+ncols*k+l)=*(parameter+ncols*k+l)*distributed_values[(int)variable_partition[k]];
                     }
                  }
               }
 
            }
            else _error_("partitioning vector should be either elements or vertex sized!");
 
         }
         else{
            /*scale all the columns by the scalar matrix:*/
            if (nrows==numberofvertices || nrows==(numberofvertices+1)){
               for(k=0;k<numberofvertices;k++){
                  if (variable_partition[k]==-1)continue;
                  else{
                     for(l=0;l<ncols;l++){
                        *(parameter+ncols*k+l)=*(parameter+ncols*k+l)*distributed_values[(int)variable_partition[k]*nt+l];
                     }
                  }
               }
            }
            else if (nrows==numberofelements || nrows==(numberofelements+1)){
               for(k=0;k<numberofelements;k++){
                  if (variable_partition[k]==-1)continue;
                  else{
                     for(l=0;l<ncols;l++){
                        *(parameter+ncols*k+l)=*(parameter+ncols*k+l)*distributed_values[(int)variable_partition[k]*nt+l];
                     }
                  }
               }
 
            }
            else _error_("partitioning vector should be either elements or vertex sized!");
         }
 
         #ifdef _DEBUG_
            PetscSynchronizedPrintf(IssmComm::GetComm(),"Parameter matrix:");
            PetscSynchronizedFlush(IssmComm::GetComm());
            for(l=0;l<ncols;l++){
               PetscSynchronizedPrintf(IssmComm::GetComm()," time %i\n",l);
               PetscSynchronizedFlush(IssmComm::GetComm());
 
               for(k=0;k<numberofvertices;k++){
                  PetscSynchronizedPrintf(IssmComm::GetComm()," node %i value %g\n",k+1,*(parameter+k*ncols+l));
                  PetscSynchronizedFlush(IssmComm::GetComm());
               }
            }
            PetscSynchronizedPrintf(IssmComm::GetComm()," descriptor: %s root %s enum: %i\n",descriptor,root,StringToEnumx(root));
            PetscSynchronizedFlush(IssmComm::GetComm());
         #endif
 
         /*Update inputs using the parameter matrix: */
         if(nrows==numberofvertices || (nrows==numberofvertices+1))
            InputUpdateFromMatrixDakotax(femmodel, parameter, nrows,ncols,StringToEnumx(root), VertexEnum);
         else
            InputUpdateFromMatrixDakotax(femmodel, parameter, nrows,ncols,StringToEnumx(root), ElementEnum);
 
         /*increment i to skip the distributed values just collected: */
         i+=npart*nt-1; //careful, the for loop will add 1.
         
         /*Free allocations: */
         xDelete<double>(parameter);
         xDelete<double>(distributed_values);
      }
      else if (strncmp(descriptor,"indexed_",8)==0){
         _error_("indexed variables not supported yet!");
      }
      else if (strncmp(descriptor,"nodal_",8)==0){
         _error_("nodal variables not supported yet!");
      }
      else{
         /*Ok, standard variable, just update inputs using the variable: */
         InputUpdateFromConstantx(femmodel,variables[i],StringToEnumx(descriptor));
      }
      variablecount++;
   }
 
   /*Free ressources:*/
   for(i=0;i<variable_partitions_num;i++){
      IssmDouble* matrix=variable_partitions[i];
      xDelete<IssmDouble>(matrix);
   }
   xDelete<IssmDouble*>(variable_partitions);
   xDelete<int>(variable_partitions_npart);
   xDelete<int>(variable_partitions_nt);
 
}