Index: /issm/trunk-jpl/src/c/cores/transient_core.cpp
===================================================================
--- /issm/trunk-jpl/src/c/cores/transient_core.cpp	(revision 25510)
+++ /issm/trunk-jpl/src/c/cores/transient_core.cpp	(revision 25511)
@@ -316,12 +316,16 @@
 	/*__________________________________________________________________________________*/
 
-	/*Get X (control) and Y (model state)*/
+	/*Get X (control)*/
 	IssmDouble *X = NULL; int Xsize;
 	GetVectorFromControlInputsx(&X,&Xsize,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,"value");
-	IssmDouble *Y = NULL; int Ysize;
-	_error_("don't know how to get Y and Ysize here...");
+
+	/*Get Y (model state) size*/
+	CountDoublesFunctor* marshallhandle1 = new CountDoublesFunctor();
+	femmodel->Marshall(marshallhandle1);
+	int Ysize = marshallhandle1->DoubleCount();
+	delete marshallhandle1;
 
 	/*Initialize Xb, Yb and Yin*/
-	double *Xb = xNewZeroInit<double>(Xsize);
+	double *Xb  = xNewZeroInit<double>(Xsize);
 	double *Yb  = xNewZeroInit<double>(Ysize);
 	int    *Yin = xNewZeroInit<int>(Ysize);
@@ -332,5 +336,7 @@
 
 	/*Reverse dependent (f)*/
-	for(int i=0; i < Ysize; i++) tape_codi.registerInput(Y[i]);
+	RegisterInputFunctor* marshallhandle2 = new RegisterInputFunctor(Yin);
+	femmodel->Marshall(marshallhandle2);
+	delete marshallhandle2;
 	for(int i=0; i < Xsize; i++) tape_codi.registerInput(X[i]);
 	SetControlInputsFromVectorx(femmodel,X);
@@ -351,5 +357,7 @@
 
 	for(int i=0;i<Xsize;i++) Xb[i] += X[i].gradient();
-	for(int i=0;i<Ysize;i++) Yb[i]  = Y[i].gradient();
+	InitAdjointFunctor* marshallhandle3 = new InitAdjointFunctor(Yb);
+	femmodel->Marshall(marshallhandle3);
+	delete marshallhandle3;
 
 	/*reverse loop for transient step (G)*/
@@ -362,8 +370,7 @@
 
 		/*We need to store the CoDiPack identifier here, since y is overwritten.*/
-		for(int i=0; i<Ysize; i++) {
-			tape_codi.registerInput(Y[i]);
-			Yin[i] = Y[i].getGradientData();
-		}
+		RegisterInputFunctor* marshallhandle4 = new RegisterInputFunctor(Yin);
+		femmodel->Marshall(marshallhandle4);
+		delete marshallhandle4;
 
 		/*Tell codipack that X is the independent*/
@@ -373,6 +380,9 @@
 		/*Get New state*/
 		transient_step(femmodel);
-		_error_("don't know how to get new Y again here...");
-		for(int i=0; i<Ysize; i++) tape_codi.registerOutput(Y[i]);
+
+		/*Register output*/
+		RegisterOutputFunctor* marshallhandle5 = new RegisterOutputFunctor();
+		femmodel->Marshall(marshallhandle5);
+		delete marshallhandle5;
 
 		/*stop tracing*/
@@ -381,9 +391,13 @@
 		/*Reverse transient step (G)*/
 		/* Using y_b here to seed the next reverse iteration there y_b is always overwritten*/
-		for(int i=0; i<Ysize; i++) Y[i].gradient() = Yb[i]; 
+		SetAdjointFunctor* marshallhandle6 = new SetAdjointFunctor(Yb);
+		femmodel->Marshall(marshallhandle6);
+		delete marshallhandle6;
+
 		tape_codi.evaluate();
+
 		/* here we access the gradient data via the stored identifiers.*/
-		for(int i=0; i<Ysize; i++) Xb[i]  = tape_codi.gradient(Yin[i]);
-		for(int i=0; i<Xsize; i++) Yb[i] += X[i].gradient();
+		for(int i=0; i<Ysize; i++) Yb[i]  = tape_codi.gradient(Yin[i]);
+		for(int i=0; i<Xsize; i++) Xb[i] += X[i].gradient();
 	}
 
@@ -398,5 +412,4 @@
 	xDelete<IssmDouble>(X);
 	xDelete<double>(Xb);
-	xDelete<IssmDouble>(Y);
 	xDelete<double>(Yb);
 	xDelete<int>(Yin);
Index: /issm/trunk-jpl/src/c/shared/io/Marshalling/Marshalling.h
===================================================================
--- /issm/trunk-jpl/src/c/shared/io/Marshalling/Marshalling.h	(revision 25510)
+++ /issm/trunk-jpl/src/c/shared/io/Marshalling/Marshalling.h	(revision 25511)
@@ -20,4 +20,5 @@
 	AD_REGISTERINPUT,
 	AD_REGISTEROUTPUT,
+	AD_INITADJOINT,
 	AD_SETADJOINT,
 #endif
@@ -152,5 +153,4 @@
 			this->tape_codi    = &(IssmDouble::getGlobalTape());
 		}
-		int DoubleCount(void){return this->double_count;};
 		void Echo(void){
 			printf("RegisterInputFunctor Echo:\n");
@@ -185,9 +185,8 @@
 
 	public:
-		RegisterOutputFunctor(void) : MarshallHandle(AD_REGISTERINPUT){
+		RegisterOutputFunctor(void) : MarshallHandle(AD_REGISTEROUTPUT){
 			this->double_count = 0;
 			this->tape_codi    = &(IssmDouble::getGlobalTape());
 		}
-		int DoubleCount(void){return this->double_count;};
 		void Echo(void){
 			printf("RegisterOutputFunctor Echo:\n");
@@ -208,4 +207,76 @@
 				for(int i=0;i<size;i++){
 					this->tape_codi->registerOutput(value[i]);
+					this->double_count++;
+				}
+			}
+		}
+}; /*}}}*/
+class InitAdjointFunctor:     public MarshallHandle{ /*{{{*/
+
+	private:
+		int                   double_count;
+		IssmDouble::TapeType* tape_codi;
+		double*               adjoint;
+
+	public:
+		InitAdjointFunctor(double* adjoint_in) : MarshallHandle(AD_INITADJOINT){
+			this->double_count = 0;
+			this->tape_codi    = &(IssmDouble::getGlobalTape());
+			this->adjoint      = adjoint_in;
+		}
+		void Echo(void){
+			printf("InitAdjointFunctor Echo:\n");
+			printf("   double_count: %i\n",double_count);
+		}
+		template<typename T> void call(T & value){
+			/*General case: do nothing*/
+		}
+		template<typename T> void call(T* & value,int size){
+			/*General case: do nothing*/
+		}
+		void call(IssmDouble value){
+			value.gradient() = adjoint[this->double_count];
+			this->double_count++;
+		}
+		void call(IssmDouble* value,int size){
+			if(value){
+				for(int i=0;i<size;i++){
+					value[i].gradient() = adjoint[this->double_count];
+					this->double_count++;
+				}
+			}
+		}
+}; /*}}}*/
+class SetAdjointFunctor:      public MarshallHandle{ /*{{{*/
+
+	private:
+		int                   double_count;
+		IssmDouble::TapeType* tape_codi;
+		double*               adjoint;
+
+	public:
+		SetAdjointFunctor(double* adjoint_in) : MarshallHandle(AD_SETADJOINT){
+			this->double_count = 0;
+			this->tape_codi    = &(IssmDouble::getGlobalTape());
+			this->adjoint      = adjoint_in;
+		}
+		void Echo(void){
+			printf("SetAdjointFunctor Echo:\n");
+			printf("   double_count: %i\n",double_count);
+		}
+		template<typename T> void call(T & value){
+			/*General case: do nothing*/
+		}
+		template<typename T> void call(T* & value,int size){
+			/*General case: do nothing*/
+		}
+		void call(IssmDouble value){
+			value.gradient() = this->adjoint[this->double_count];
+			this->double_count++;
+		}
+		void call(IssmDouble* value,int size){
+			if(value){
+				for(int i=0;i<size;i++){
+					value[i].gradient() = this->adjoint[this->double_count];
 					this->double_count++;
 				}
@@ -221,8 +292,9 @@
 		case MARSHALLING_SIZE: {SizeCheckpointFunctor*  temp = xDynamicCast<SizeCheckpointFunctor*>(this);  temp->call(value); break;}
 #ifdef _HAVE_CODIPACK_
-		case AD_COUNTDOUBLES:  {CountDoublesFunctor*    temp = xDynamicCast<CountDoublesFunctor*>(this);    temp->call(value); break;}
-		case AD_REGISTERINPUT: {RegisterInputFunctor*   temp = xDynamicCast<RegisterInputFunctor*>(this);   temp->call(value); break;}
-		case AD_REGISTEROUTPUT:{RegisterOutputFunctor*  temp = xDynamicCast<RegisterOutputFunctor*>(this);  temp->call(value); break;}
-		//case AD_SETADJOINT:    {SetAdjointFunction*     temp = xDynamicCast<SetAdjointFunction*>(this);     temp->call(value); break;}
+		case AD_COUNTDOUBLES:  {CountDoublesFunctor*   temp = xDynamicCast<CountDoublesFunctor*>(this);    temp->call(value); break;}
+		case AD_REGISTERINPUT: {RegisterInputFunctor*  temp = xDynamicCast<RegisterInputFunctor*>(this);   temp->call(value); break;}
+		case AD_REGISTEROUTPUT:{RegisterOutputFunctor* temp = xDynamicCast<RegisterOutputFunctor*>(this);  temp->call(value); break;}
+		case AD_INITADJOINT:   {InitAdjointFunctor*    temp = xDynamicCast<InitAdjointFunctor*>(this);     temp->call(value); break;}
+		case AD_SETADJOINT:    {SetAdjointFunctor*     temp = xDynamicCast<SetAdjointFunctor*>(this);      temp->call(value); break;}
 #endif
 		default: _error_("Operation "<<OperationNumber()<<" not supported yet");
@@ -235,8 +307,9 @@
 		case MARSHALLING_SIZE: {SizeCheckpointFunctor*  temp = xDynamicCast<SizeCheckpointFunctor*>(this);  temp->call(value,size); break;}
 #ifdef _HAVE_CODIPACK_
-		case AD_COUNTDOUBLES:  {CountDoublesFunctor*    temp = xDynamicCast<CountDoublesFunctor*>(this);    temp->call(value,size); break;}
-		case AD_REGISTERINPUT: {RegisterInputFunctor*   temp = xDynamicCast<RegisterInputFunctor*>(this);   temp->call(value,size); break;}
-		case AD_REGISTEROUTPUT:{RegisterOutputFunctor*  temp = xDynamicCast<RegisterOutputFunctor*>(this);  temp->call(value,size); break;}
-		//case AD_SETADJOINT:    {SetAdjointFunction*     temp = xDynamicCast<SetAdjointFunction*>(this);     temp->call(value,size); break;}
+		case AD_COUNTDOUBLES:  {CountDoublesFunctor*   temp = xDynamicCast<CountDoublesFunctor*>(this);    temp->call(value,size); break;}
+		case AD_REGISTERINPUT: {RegisterInputFunctor*  temp = xDynamicCast<RegisterInputFunctor*>(this);   temp->call(value,size); break;}
+		case AD_REGISTEROUTPUT:{RegisterOutputFunctor* temp = xDynamicCast<RegisterOutputFunctor*>(this);  temp->call(value,size); break;}
+		case AD_INITADJOINT:   {InitAdjointFunctor*    temp = xDynamicCast<InitAdjointFunctor*>(this);     temp->call(value,size); break;}
+		case AD_SETADJOINT:    {SetAdjointFunctor*     temp = xDynamicCast<SetAdjointFunctor*>(this);      temp->call(value,size); break;}
 #endif
 		default: _error_("Operation "<<OperationNumber() <<" not supported yet");