ivp_buoyancy_solver.cxx

hl2 source code. Do not use it illegal.

	IVP_U_Float_Point area_center_total;	        //area center of the whole triangle
	area_center_total.add(triangle_points[0], triangle_points[1]);
	area_center_total.add(triangle_points[2]);
		
	//compute the area center of the cut-off triangle above the surface (Ce2 = (SP1+SP2+C)/3)
	IVP_U_Float_Point area_center_part_above;	//area center of the triangle's part above the surface
	area_center_part_above.set_multiple(triangle_points[index_neg], 3.0f - factor1-factor2);
	area_center_part_above.add_multiple(triangle_points[index_neg+1], factor1);
	area_center_part_above.add_multiple(triangle_points[index_neg+2], factor2);
		
		
	//compute the area center of the cut-off quadrangle under the surface
	//Ce1 = (C - Ce2 * k) / (1-k)
	IVP_DOUBLE k = factor1 * factor2;		//factor the area center is moved
	IVP_U_Float_Point area_center_part_under;	//area center of the triangle's part under the surface
	area_center_part_above.mult(k);
	area_center_part_under.set(&area_center_total);
	area_center_part_under.subtract(&area_center_part_above);
	if ((1.0f-k) < buoyancy_eps) {  //both other vertices are exactly on the surface, so there is no volume under the surface
	    volume_under = 0.0f;
	    volume_center_under.set(0.0f, 0.0f, 0.0f);
	    return;
	}
	area_center_part_under.mult( 1.0f / (1.0f-k));
		
	//compute the volume center of the part under the surface from area_center_part_under
	area_center_part_under.add(s_point);
		
	area_center_part_under.mult(volume_pyramid - volume_pyramid_above);
	volume_center_under.add(&area_center_part_under);
    }
    break;
    case 2: { //two vertices are above the surface, the vertex with index "index_positiv" is under the surface
		
		
	//compute the first disection point on the surface
	IVP_DOUBLE factor1 = distance[index_positiv] / (distance[index_positiv] - distance[index_positiv+1]);  //division by zero not possible, because distance[index_positiv+1] < 0 and buoyancy_eps is added to each
		
	//compute the second disection point on the surface
	IVP_DOUBLE factor2 = distance[index_positiv] / (distance[index_positiv] - distance[index_positiv+2]);  //division by zero not possible, because distance[index_positiv+2] < 0 and buoyancy_eps is added to each
		
		
	IVP_IF(1){
	    //for debugging reasons - begin
	    IVP_U_Float_Point sp1, sp2;
	    IVP_U_Point sp1_ws, sp2_ws;
	    sp1.set_interpolate(triangle_points[index_positiv], triangle_points[index_positiv+1], factor1);
	    sp2.set_interpolate(triangle_points[index_positiv], triangle_points[index_positiv+2], factor2);
	    IVP_Cache_Object *cache_object = object->get_cache_object();  //needed to transform coordinates into other coord. systems
	    cache_object->transform_position_to_world_coords(&sp1, &sp1_ws);
	    cache_object->transform_position_to_world_coords(&sp2, &sp2_ws);
	    cache_object->remove_reference();
			
	    IVP_U_Float_Point edge12;
	    edge12.subtract(&sp1_ws, &sp2_ws);
	    IVP_U_Point tmp_sp2_ws;
	    tmp_sp2_ws.set(&sp2_ws);
	    environment->add_draw_vector(&tmp_sp2_ws, &edge12,"",4);
	    //for debugging reasons - end
	}
		
	//compute volume of whole pyramid triangle is part of
	IVP_DOUBLE volume_pyramid = compute_pyramid_volume( triangle_points[0], triangle_points[1], triangle_points[2], s_point );
		
	//compute volume of pyramid under the surface
	IVP_DOUBLE volume_pyramid_under = factor1 * factor2 * volume_pyramid;
		
	volume_under += volume_pyramid_under;
		
		
	/*** compute the area centers ***/
		
	//compute the area center of the cut-off triangle under the surface (Ce2 = (SP1+SP2+C)/3)
	IVP_U_Float_Point area_center_part_under;
	area_center_part_under.set_multiple(triangle_points[index_positiv], 3.0f-factor1-factor2);
	area_center_part_under.add_multiple(triangle_points[index_positiv+1], factor1);
	area_center_part_under.add_multiple(triangle_points[index_positiv+2], factor2);
		
	//compute the volume center of the pyramid under the surface from area_center_part_under
	area_center_part_under.add(s_point);
		
	//save results to class variables
	area_center_part_under.mult(volume_pyramid_under); //weight the volume center of this part of the object with its volume
	volume_center_under.add(&area_center_part_under);
    }
    break;
    default:
	break;
    }
	
}


void IVP_Buoyancy_Solver::compute_rotation_and_translation_values_for_one_triangle(IVP_Real_Object *object,
										   const IVP_Compact_Triangle *current_triangle,
										   const IVP_U_Float_Point *triangle_points[5],
										   const IVP_Compact_Ledge *current_ledge,
										   const IVP_FLOAT distance[6],
										   const int &decision,
										   const int &index_positiv,
										   const int &index_neg) {
	
	
	
    //declare some common variables
    IVP_U_Float_Point area_center_part_under_os;
    IVP_FLOAT triangle_surface_under_factor;
    
    switch (decision) {
    case 0: { //current_triangle is completely under the surface
		
	//compute the area center of the current_triangle
	area_center_part_under_os.add(triangle_points[0], triangle_points[1]);
	area_center_part_under_os.add(triangle_points[2]);
	area_center_part_under_os.mult( 1.0f/3.0f );
		
	triangle_surface_under_factor = 1.0f;
		
	break; 
    }
		
    case 1: { //the vertex with index "index_neg" is above the surface
		
	//compute the relation factor between the vertex above the surface and the first one under the surface
	IVP_DOUBLE factor1 = -distance[index_neg] / (distance[index_neg+1] - distance[index_neg]);  //division by zero not possible, because distance[index_neg] < 0 and buoyancy_eps is added to each
		
	//compute the relation factor between the vertex above the surface and the second one under the surface
	IVP_DOUBLE factor2 = -distance[index_neg] / (distance[index_neg+2] - distance[index_neg]);  //division by zero not possible, because distance[index_neg] < 0 and buoyancy_eps is added to each
		
	//compute the area center of the current_triangle
	IVP_U_Float_Point area_center_of_triangle_os;
	area_center_of_triangle_os.add(triangle_points[0], triangle_points[1]);
	area_center_of_triangle_os.add(triangle_points[2]);
	area_center_of_triangle_os.mult( 1.0f/3.0f );
		
	//compute the area center of the cut-off triangle above the surface (Ce2 = (SP1+SP2+C)/3)
	IVP_U_Float_Point area_center_part_above_os;       //area center of the triangle's part above the surface
	area_center_part_above_os.set_multiple(triangle_points[index_neg], (3.0f-factor1-factor2) * (1.0f/3.0f));
	area_center_part_above_os.add_multiple(triangle_points[index_neg+1], factor1 * (1.0f/3.0f));
	area_center_part_above_os.add_multiple(triangle_points[index_neg+2], factor2 * (1.0f/3.0f));
		
	//compute the area center of the cut-off quadrangle under the surface
	//Ce1 = (C - Ce2 * k) / (1-k)
	IVP_DOUBLE k = factor1 * factor2;           //factor the area center is moved
		
	area_center_part_above_os.mult(k);
	area_center_part_under_os.set(&area_center_of_triangle_os);
	area_center_part_under_os.subtract(&area_center_part_above_os);
	
	area_center_part_under_os.mult( 1.0f / (1.0f-k));
	if ((1.0f-k) < buoyancy_eps) {  //both other vertices are exactly on the surface, so there is nothing to damp
	    return;
	}
	triangle_surface_under_factor = 1.0f - k;
		
	break;
    }
		
    case 2: { //two vertices are above the surface, the vertex with index "index_positiv" is under the surface
		
	//compute the relation factor between the vertex under the surface and the first one above the surface
	IVP_DOUBLE factor1 = distance[index_positiv] / (distance[index_positiv] - distance[index_positiv+1]);  //division by zero not possible, because distance[index_positiv+1] < 0 and buoyancy_eps is added to each
		
	//compute the relation factor between the vertex under the surface and the first one above the surface
	IVP_DOUBLE factor2 = distance[index_positiv] / (distance[index_positiv] - distance[index_positiv+2]);  //division by zero not possible, because distance[index_positiv+2] < 0 and buoyancy_eps is added to each
		
	//compute the area center of the cut-off triangle under the surface (Ce2 = (SP1+SP2+C)/3)
	area_center_part_under_os.set_multiple(triangle_points[index_positiv], (3.0f-factor1-factor2) * (1.0f/ 3.0f));
	area_center_part_under_os.add_multiple(triangle_points[index_positiv+1], factor1 * (1.0f/3.0f));
	area_center_part_under_os.add_multiple(triangle_points[index_positiv+2], factor2 * (1.0f/3.0f));
		
	triangle_surface_under_factor = factor1*factor2;
		
	break;
    }
		
    default:
	return;
    }
    
	
    {
	// compute normal_os of triangle
	const IVP_Compact_Edge *edge = current_triangle->get_first_edge();
	IVP_U_Float_Point normal_os;
	IVP_Compact_Ledge_Solver::calc_hesse_vec_object_not_normized(edge, current_ledge, &normal_os);
		
	if (normal_os.quad_length() < buoyancy_eps) {  //savety check
	    //if this is the case then it is better not to damp this triangle
	    return;
	}
		
	//get surface speed of the object at the area center of the part under the surface
	IVP_U_Float_Point speed_of_object_os;
	ivp_core_get_surface_speed_os(core,object, &area_center_part_under_os, &speed_of_object_os);
		
	//calc the relative speed of the medium
	IVP_U_Float_Point rel_speed_of_current_os;
	rel_speed_of_current_os.subtract( &resulting_speed_of_current_os, &speed_of_object_os);
		
	//get quad_length of relative speed, needed later
	IVP_FLOAT quad_length_of_speed = rel_speed_of_current_os.quad_length();
	if (quad_length_of_speed < buoyancy_eps) {  //savety check
	    //better no dampening of this triangle
	    quad_length_of_speed = 0.0f;
	    rel_speed_of_current_os.set(0.0f, 1.0f, 0.0f);
	}

		
	IVP_FLOAT surface_factor_check = normal_os.dot_product( &rel_speed_of_current_os );
	//cancel calculation of triangle if turned away from the current
	if (surface_factor_check > buoyancy_eps) {
	    return;
	}
	
	rel_speed_of_current_os.normize();

	IVP_U_Float_Point normized_normal_os;
	normized_normal_os.set(&normal_os);
	IVP_DOUBLE len_of_unscaled_normal = normized_normal_os.real_length_plus_normize();
	IVP_FLOAT surface_factor = normized_normal_os.dot_product( &rel_speed_of_current_os );
		
	//calc the surface content of the whole triangle
	IVP_FLOAT triangle_surface_content = 0.5f * len_of_unscaled_normal;

	//calculate the projected surface content
	IVP_FLOAT visible_surface_content_under = -surface_factor * triangle_surface_under_factor * triangle_surface_content;
 
		
	IVP_DOUBLE impulse;
	//impulse = (visible_surface_content_under * quad_length_of_speed * medium_damp_factor);

	//calculate the pressure dampening
	impulse = (visible_surface_content_under * quad_length_of_speed * pressure_damp_factor);
	//calculate the friction dampening
	impulse += (triangle_surface_content * quad_length_of_speed * friction_damp_factor);
	IVP_U_Float_Point impulse_vector;

	if (simulate_wing_behavior) {
	    impulse_vector.set_multiple( &normized_normal_os, -impulse );
	} else {
	    impulse_vector.set_multiple( &rel_speed_of_current_os, impulse );
	}
		
	IVP_IF(0) {
	    IVP_U_Point center_of_triangle_ws;
	    IVP_U_Float_Point impulse_vector_ws;
	    {
		IVP_Cache_Object *cache_object = object->get_cache_object();
		cache_object->transform_position_to_world_coords(&area_center_part_under_os, &center_of_triangle_ws);
		cache_object->transform_vector_to_world_coords(&impulse_vector, &impulse_vector_ws);
		cache_object->remove_reference();
	    }
	    environment->add_draw_vector(&center_of_triangle_ws, &impulse_vector_ws,"",2);
	}
		
	//add surface content, which is under the water, to the corresponding global variable
	object_visible_surface_content_under += IVP_Inline_Math::fabsd(visible_surface_content_under);
		
	//calc cross product between impulse vector and area center of triangle and add it to the global variable
	IVP_U_Float_Point impulse_x_point;
	impulse_x_point.inline_calc_cross_product(&area_center_part_under_os, &impulse_vector);
	sum_impulse_x_point.add(&impulse_x_point);
		
		
	//calc the move-direction vector
	IVP_U_Float_Point move_direction_os;
	move_direction_os.set_orthogonal_part(&rel_speed_of_current_os, &normized_normal_os);
		
		//calc cross product between impulse_vector and move-direction vector and add it to the global variable
	IVP_U_Float_Point impulse_x_movevector;
	impulse_x_movevector.inline_calc_cross_product(&move_direction_os, &impulse_vector);
	sum_impulse_x_movevector.add(&impulse_x_movevector);
		
	//add impulse to global variable
	sum_impulse.add(&impulse_vector);
    }
	
    
	
}


void IVP_Buoyancy_Solver::compute_values_for_one_ledge(IVP_Real_Object *object,
						       const IVP_Compact_Ledge *current_ledge,
						       const IVP_U_Float_Hesse *surface_os,
						       const IVP_U_Float_Point *s_point) {
	
    int nr_of_triangles = current_ledge->get_n_triangles();  //fetch number of triangles in current_ledge
    const IVP_Compact_Triangle *triangle = current_ledge->get_first_triangle();  //fetch first triangle in current_ledge
	
    for (int i=nr_of_triangles; --i >= 0; triangle = triangle->get_next_tri()) {  //fetch next triangle of ledge
	compute_values_for_one_triangle(object, triangle, surface_os, s_point, current_ledge);
    }
}