 {
 	this->partitions = partitions;
 	return 0;
-}
-
-vector<float3> addMassMatrix(double rho, double A, double l)
-{
-	vector < float3 > massNew;
-	float3 massEntry;
-	massEntry.x = 0;
-	massEntry.y = 0;
-	massEntry.z = 0.13e2 / 0.35e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 0;
-	massEntry.y = 3;
-	massEntry.z = 0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 0;
-	massEntry.y = 6;
-	massEntry.z = 0.9e1 / 0.70e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 0;
-	massEntry.y = 9;
-	massEntry.z = -0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 1;
-	massEntry.y = 1;
-	massEntry.z = 0.13e2 / 0.35e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 1;
-	massEntry.y = 4;
-	massEntry.z = 0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 1;
-	massEntry.y = 7;
-	massEntry.z = 0.9e1 / 0.70e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 1;
-	massEntry.y = 10;
-	massEntry.z = -0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 2;
-	massEntry.y = 2;
-	massEntry.z = 0.13e2 / 0.35e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 2;
-	massEntry.y = 5;
-	massEntry.z = 0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 2;
-	massEntry.y = 8;
-	massEntry.z = 0.9e1 / 0.70e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 2;
-	massEntry.y = 11;
-	massEntry.z = -0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 3;
-	massEntry.y = 0;
-	massEntry.z = 0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 3;
-	massEntry.y = 3;
-	massEntry.z = rho * A * l * l * l / 0.105e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 3;
-	massEntry.y = 6;
-	massEntry.z = 0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 3;
-	massEntry.y = 9;
-	massEntry.z = -rho * A * l * l * l / 0.140e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 4;
-	massEntry.y = 1;
-	massEntry.z = 0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 4;
-	massEntry.y = 4;
-	massEntry.z = rho * A * l * l * l / 0.105e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 4;
-	massEntry.y = 7;
-	massEntry.z = 0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 4;
-	massEntry.y = 10;
-	massEntry.z = -rho * A * l * l * l / 0.140e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 5;
-	massEntry.y = 2;
-	massEntry.z = 0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 5;
-	massEntry.y = 5;
-	massEntry.z = rho * A * l * l * l / 0.105e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 5;
-	massEntry.y = 8;
-	massEntry.z = 0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 5;
-	massEntry.y = 11;
-	massEntry.z = -rho * A * l * l * l / 0.140e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 6;
-	massEntry.y = 0;
-	massEntry.z = 0.9e1 / 0.70e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 6;
-	massEntry.y = 3;
-	massEntry.z = 0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 6;
-	massEntry.y = 6;
-	massEntry.z = 0.13e2 / 0.35e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 6;
-	massEntry.y = 9;
-	massEntry.z = -0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 7;
-	massEntry.y = 1;
-	massEntry.z = 0.9e1 / 0.70e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 7;
-	massEntry.y = 4;
-	massEntry.z = 0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 7;
-	massEntry.y = 7;
-	massEntry.z = 0.13e2 / 0.35e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 7;
-	massEntry.y = 10;
-	massEntry.z = -0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 8;
-	massEntry.y = 2;
-	massEntry.z = 0.9e1 / 0.70e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 8;
-	massEntry.y = 5;
-	massEntry.z = 0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 8;
-	massEntry.y = 8;
-	massEntry.z = 0.13e2 / 0.35e2 * rho * A * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 8;
-	massEntry.y = 11;
-	massEntry.z = -0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 9;
-	massEntry.y = 0;
-	massEntry.z = -0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 9;
-	massEntry.y = 3;
-	massEntry.z = -rho * A * l * l * l / 0.140e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 9;
-	massEntry.y = 6;
-	massEntry.z = -0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 9;
-	massEntry.y = 9;
-	massEntry.z = rho * A * l * l * l / 0.105e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 10;
-	massEntry.y = 1;
-	massEntry.z = -0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 10;
-	massEntry.y = 4;
-	massEntry.z = -rho * A * l * l * l / 0.140e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 10;
-	massEntry.y = 7;
-	massEntry.z = -0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 10;
-	massEntry.y = 10;
-	massEntry.z = rho * A * l * l * l / 0.105e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 11;
-	massEntry.y = 2;
-	massEntry.z = -0.13e2 / 0.420e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 11;
-	massEntry.y = 5;
-	massEntry.z = -rho * A * l * l * l / 0.140e3;
-	massNew.push_back(massEntry);
-	massEntry.x = 11;
-	massEntry.y = 8;
-	massEntry.z = -0.11e2 / 0.210e3 * rho * A * l * l;
-	massNew.push_back(massEntry);
-	massEntry.x = 11;
-	massEntry.y = 11;
-	massEntry.z = rho * A * l * l * l / 0.105e3;
-	massNew.push_back(massEntry);
-
-	return massNew;
 }
 
 int ANCFSystem::addParticle(Particle* particle)

	fext_h.push_back(-rho * A * a * a * GRAVITYx / 0.12e2);
	fext_h.push_back(-rho * A * a * a * GRAVITYy / 0.12e2);
	fext_h.push_back(-rho * A * a * a * GRAVITYz / 0.12e2);
	
//	//add mass matrix
//	vector<float3> mass = addMassMatrix(rho,A,a);
//
//	for(int i=0;i<mass.size();i++)
//	{
//		massI_h.push_back(mass[i].x+12*(elements.size()-1));
//		massJ_h.push_back(mass[i].y+12*(elements.size()-1));
//		mass_h.push_back(static_cast<double>(mass[i].z));
//	}

	for(int i=0;i<12;i++)
	{

	return 0;
}

int ANCFSystem::initializeSystem()
{
	ANCFSystem::updatePhiq();
	ANCFSystem::calculateInitialPhi();

	for(int i=0;i<constraints.size();i++)
	{
