finallogdet.cpp

function_that_estimates_logdet
{
	CSparseMatrixOperator<float64_t> *op
	CSerialComputationEngine *e
	CLogRationalApproximationCOCGM log_op(op, e) 
	// or CLogRationalApproximationIndividual log_op(op, e)
	CProbingSampler *sampler
	CLogDetEstimator est(log_op, e, sampler)
	SGVector<float64_t> samples=est.sample(num_logdet_estimates)
}

CLogDetEstimator::sample(index_t num_logdet_estimates)
{
	//calls eigensolver, computes shift, weight and integration const
	m_operator->init()
	CDynamicObjectArray job_result_aggregations;
	for i=0 to num_logdet_estimates
	{
		//returns the number of colors in the graph
		N=m_trace_sampler->get_num_samples() 
		// assuming that Dan says yes to our approach
		for j=0 to N
		{
			SGVector<float64_t> s=m_trace_sampler->sample(j)
			CJobResultAggregator* agg=m_operator_log->submit_jobs(s)
			// we always get one aggregator per vector, 
			// no matter how many jobs are there per vector
			job_result_aggregations.appened_element(agg)
		}
	}

	m_independent_computaion_engine->wait_for_all()

	for all elements in job_result_aggregations
	{
		CJobResultAggregator *agg=job_result_aggregations.get_element(i)
		agg->finalize()
		float64_t result=agg->get_result()
		// do something with the result, running avg, avg etc
	}
}

CLogRationalApproximation<float64_t>::init()
{
	// compute shifts, weights, const, common for all
}

CLogRationalApproximationCOCGM::submit_jobs(SGVector<float64_t> s)
{
	CStoreScalarAggregator<float64_t> *agg
	CLogRationalApproximationCOCGMJob *job(agg) // just one job per vector
	//hands over all responsibilities of the job to the engine
	m_computation_engine->submit_job(job)
	return agg
}

CLogRationalApproximationIndividual::submit_jobs(SGVector<float64_t> s)
{
	CIndividualJobResultAggregator *agg(s)
	for i=0 to num_shifts
	{
		CLogRationalApproximationIndividualJob *job(agg, i) 
		// #shifts jobs per vector, same aggregator, i=index
		m_computation_engine->submit_job(job)
	}
	return agg
}

CSerialComputationEngine::submit_job(CIndependentJob *job)
{
	job->compute() //solve systems, therefore blocks until compute returns
}

CSerialComputationEngine::wait_for_all()
{
	//does nothing
}

CLogRationalApproximationCOCGMJob::compute()
{
	// use the solver in CLogRationalApproximationCOCGM to solve system
	SGVector<complex64_t> vec=CLogRationalApproximation::m_linear_solver->solve(m_operator, m_vector)
	// don't have to pass weights, shfits, const and the operator
	m_job_aggregator->submit_result(CScalarResult<float64_t>(m_vector.dot(imag(vec))))
}

CLogRationalApproximationIndividualJob::compute()
{
	//use the index to form this, using shifts, weights etc
	CSpasrseLinearOperator<complex64_t> shifted_op; 
	SGVector<complex64_t> vec=CLogRationalApproximation::m_linear_solver->solve(shifted_op, m_vector)
	m_job_aggregator->submit_result(CVectorResult<complex64_t>(vec))
}

CIndividualJobResultAggregator::submit_result(CVectorResult<complex64_t> v)
{
	// adds the param v (newcomer solution vector) to the already existing sum
	// which will then be used in the vec-vec product in the finalize
	// will be called #shifts times for individual solve
	m_aggregate+=v.get_result() // all complex
}

CIndividualJobResultAggregator::finalize()
{
	//s is stored inside the aggregator, increased refcount
	float64_t result=s.dot(get_imag(m_aggregate))
	// store the scalar result inside it, get_aggregation returns it
	m_result=CScalarResult<float64_t>(result)
	set_some_flag_true()
}

CStoreScalarAggregator<float64_t>::submit_result(CScalarResult<float64_t> v)
{
	m_aggregate+=v.get_result()
	// although not needed for COCG_M, will be called only once with each vector
	// so its just once, but written like this for general purpose
}

CStoreScalarAggregator<float6464_t>::finalize()
{
	m_result=CScalarResult<float64_t>(m_aggregate)
	set_some_flag_true() // does nothing basiclly
}

CJobResultAggregator::get_aggregation()
{
	return m_result->get_result()
}