tcstan2.cpp

并行TIN生成算法, 基于DeWall算法理论实现

			return;
		}
	}
}

// 批量单元搜索
TCS_V00 FndEmp2(TIN_G *G, TCS_F64 *P, TCS_F64 *p1, TCS_F64 *p2, TCS_F64 *ce, TCS_F64 &d2, CTcsBL *L, TCS_S32 co, TCS_S32 *c, TCS_S32 *r)
{
	TCS_S32 i, k, ci;

	if (r[2] != r[0])	// 即r[2] != r[4]
	{
		ci = r[0] * co;
		for (k = c[0] + 1; k < c[1]; k++)
		{
			FndEmp2(G, P, p1, p2, ci + k, ce, d2, L);
			if (L->GetSize() > 0)
				return;
		}
	}
	if (r[3] != r[1])	// 即r[3] != r[5]
	{
		ci = r[1] * co;
		for (k = c[0] + 1; k < c[1]; k++)
		{
			FndEmp2(G, P, p1, p2, ci + k, ce, d2, L);
			if (L->GetSize() > 0)
				return;
		}
	}
	if (c[2] != c[0])	// 即c[2] != c[4]
	{
		for (i = r[0] + 1; i < r[1]; i++)
		{
			ci = i * co + c[0];
			FndEmp2(G, P, p1, p2, ci, ce, d2, L);
			if (L->GetSize() > 0)
				return;
		}
	}
	if (c[3] != c[1])	// 即c[3] != c[5]
	{
		for (i = r[0] + 1; i < r[1]; i++)
		{
			ci = i * co + c[1];
			FndEmp2(G, P, p1, p2, ci, ce, d2, L);
			if (L->GetSize() > 0)
				return;
		}
	}
	// 处理4个角单元
	if (c[2] != c[4] || r[2] != r[4])
	{
		ci = r[0] * co + c[0];
		FndEmp2(G, P, p1, p2, ci, ce, d2, L);
		if (L->GetSize() > 0)
			return;
	}
	if (c[2] != c[4] || r[3] != r[5])
	{
		ci = r[1] * co + c[0];
		FndEmp2(G, P, p1, p2, ci, ce, d2, L);
		if (L->GetSize() > 0)
			return;
	}
	if (c[3] != c[5] || r[2] != r[4])
	{
		ci = r[0] * co + c[1];
		FndEmp2(G, P, p1, p2, ci, ce, d2, L);
		if (L->GetSize() > 0)
			return;
	}
	if (c[3] != c[5] || r[3] != r[5])
	{
		ci = r[1] * co + c[1];
		FndEmp2(G, P, p1, p2, ci, ce, d2, L);
		if (L->GetSize() > 0)
			return;
	}
}


// 生成构网分区
//#define TEST_SIN
TCS_V00 CnsTen2(TAN_PAR *TAN)
{
	TIN_G *G = TAN->G;
	TCS_F64 *P = TAN->P;
	TCS_S32 *vpn = TAN->vpn;
	TCS_S32 *nzn = TAN->nzn;

	TCS_F64 d2, ce[2];
	ce[0] = (G->bnd[0] + G->bnd[2]) / 2;
	ce[1] = (G->bnd[1] + G->bnd[3]) / 2;

	TCS_S32 sid[2], eid[2], sip[2], eip[2];
	// 搜索离区域中心最近的点
	FndNep2(G, P, ce, sid[0], sid[1]);
	TCS_F64 *sp = P + 2 * sid[0];	sip[0] = sid[0], sip[1] = sid[1];
	// 搜索第2点:最近法则
	FndNep2(G, P, sp, eid[0], eid[1]);
	TCS_F64 *ep = P + 2 * eid[0];	eip[0] = eid[0], eip[1] = eid[1];
	// 搜索第3点:空圆法则
	CTcsBL *ARL = TAN->ARL;	// 隔离区链表
	CTcsBL *ATL = TAN->ATL;	// 三角形链表
	CTcsHL *BSL = TAN->BSL;	// 静止边链表
	CTcsBL *APL = new CTcsBL;	// 活动点链表
	CTcsSL *SSL = new CTcsSL;	// 候选边链表
	CTcsHL *ASL = new CTcsHL(128, SIDE_EQLH, SIDE_HASH);	// 活动边链表
	{
		ce[0] = (sp[0] + ep[0]) / 2;
		ce[1] = (sp[1] + ep[1]) / 2;
		// 计算初始搜索范围
		TCS_S32 i, k, c[6], r[6];
		TCS_F64 d = sqrt(EuLen2D(sp, ce));
		ComBnd2(G, ce, d, c[0], c[1], r[0], r[1]);

		TCS_S32 ci, co = G->num[0], ro = G->num[1];

		for (i = r[0]; i <= r[1]; i++)
		{
			ci = i * co;
			for (k = c[0]; k <= c[1]; k++)
			{
				FndEmp2(G, P, sp, ep, ci + k, ce, d2, APL);
				if (APL->GetSize() > 0)
					break;
			}
			if (APL->GetSize() > 0)
				break;
		}

		if (APL->GetSize() == 0)
		{
			// 继续搜索
			c[2] = c[0];	c[3] = c[1];		r[2] = r[0];	r[3] = r[1];
			c[4] = 0;		c[5] = co - 1;		r[4] = 0;		r[5] = ro - 1;
			while (1)
			{
				c[0]--;	if (c[0] < c[4])	c[0] = c[4];
				c[1]++;	if (c[1] > c[5])	c[1] = c[5];
				r[0]--;	if (r[0] < r[4])	r[0] = r[4];
				r[1]++;	if (r[1] > r[5])	r[1] = r[5];
				FndEmp2(G, P, sp, ep, ce, d2, APL, co, c, r);
				if (APL->GetSize() > 0)
					break;	// 已经找到

				c[2] = c[0], c[3] = c[1], r[2] = r[0], r[3] = r[1];
				if (c[2] == c[4] && c[3] == c[5] &&
					r[2] == r[4] && r[3] == r[5])
					break;
			}
		}
		if (APL->GetSize() == 0)	// 终止
		{
			delete ASL;	delete APL;
			delete SSL;
			return;
		}
		// 共圆点弧位排序
		if (APL->GetSize() > 1)
			ArcPos2(APL, P, ce, d2);
		// 逐点构三角形, 如4点以上共圆需局部优化
		ConTen2(G, P, sid, eid, APL, ATL, ASL, BSL, SSL);
	}
	// 单核测试
#ifdef TEST_SIN
	TIN_S *S = (TIN_S*)ASL->GetHead();
	while (S)
	{
		FndEmp2(G, S, P, ce, d2, APL);
		if (APL->GetSize() == 0)
		{
			// 减少计数器
			G->pid[S->sid[0]]--;	G->pid[S->eid[0]]--;
			// 删除归零点
			if (G->pid[S->sid[0]] == 0)	{	RmvPid(G->idx + S->sid[1], S->sid[0]);	G->vpn++;	}
			if (G->pid[S->eid[0]] == 0)	{	RmvPid(G->idx + S->eid[1], S->eid[0]);	G->vpn++;	}
			// 收集边对象
			S->pid[1] = -1;		ASL->DetData(S);	BSL->PutData(S, S);
		}
		else
		{
			// 共圆点排序
			if (APL->GetSize() > 1)
				ArcPos2(APL, P, ce, d2);
			// 构造三角形
			sid[0] = S->sid[0], sid[1] = S->sid[1];
			eid[0] = S->eid[0], eid[1] = S->eid[1];
			ConTen2(G, P, sid, eid, APL, ATL, ASL, BSL);
		}

		S = (TIN_S*)ASL->GetHead();
	}
	// 释放内存资源
	delete APL;	delete ASL;	delete SSL;
	return;
#endif
	// 生成构网分区
	ce[0] = (sp[0] + ep[0]) / 2;
	ce[1] = (sp[1] + ep[1]) / 2;
	if (sp[0] == ep[0])
	{
		// 从Y分隔线开始
		TIN_S *S = (TIN_S*)SSL->PopData();
		while (S)
		{
			// 如果是活动边, 判断是否与Y分隔线相交
			if (ASL->GetData(S) && LySide2(ce[1], P, S))
				ConTyn2(G, S, P, ce[1], ATL, ASL, BSL);
			S = (TIN_S*)SSL->PopData();
		}
		// 构造隔离区
		TIN_R *UR = new TIN_R;	// 上隔离区
		TIN_R *DR = new TIN_R;	// 下隔离区

		ARL->AddTail(UR);	ARL->AddTail(DR);
		ConRyn2(G, UR, DR, ASL, P, ce[0], ce[1]);
	}
	else
	{
		// 从X分隔线开始
		TIN_S *S = (TIN_S*)SSL->PopData();
		while (S)
		{
			// 如果是活动边, 判断是否与X分隔线相交
			if (ASL->GetData(S) && LxSide2(ce[0], P, S))
				ConTxn2(G, S, P, ce[0], ATL, ASL, BSL);
			S = (TIN_S*)SSL->PopData();
		}
		// 构造隔离区
		TIN_R *LR = new TIN_R;	// 左隔离区
		TIN_R *RR = new TIN_R;	// 右隔离区

		ARL->AddTail(LR);	ARL->AddTail(RR);
		ConRxn2(G, LR, RR, ASL, P, ce[0], ce[1]);
	}
	// 释放内存资源
	delete APL;	delete ASL;	delete SSL;

	// 按点数进行循环分区
	// 按点数进行分区
	TCS_S32 zn = *vpn / *nzn;
	if (zn == 0)	zn = 1;
	// 采用2分法, 即1分2
	CTcsBL *B1 = ARL;
	CTcsBL *B2 = new CTcsBL;

	// 调度线程
	TCS_V00 **nven = new TCS_V00*[TAN->unnu];
	TUN_PAR *unpa = (TUN_PAR*)TAN->unpa;
	for (long i = 0; i < TAN->unnu; i++)
		nven[i] = unpa[i].nven[1];
	// 循环构造
	while (B1->GetSize() < zn)
	{
		TIN_R *R = (TIN_R*)B1->RemoveHead();
		while (R)
		{
			// 等待线程空闲
			TCS_U32 idx = ::WaitForMultipleObjects(TAN->unnu, nven, 0, INFINITE);
			if (idx >= WAIT_OBJECT_0 && idx < WAIT_OBJECT_0 + TAN->unnu)
			{
				idx = idx - WAIT_OBJECT_0;
				unpa[idx].P = P;
				unpa[idx].R = R;
				::SetEvent(unpa[idx].nven[0]);
			}
			B2->AddTail(R);
			// 如果中途终止
			if (*(TAN->exit) == 1)
				break;

			R = (TIN_R*)B1->RemoveHead();
		}
		// 等待本次隔离任务完成
		for (long i = 0; i < TAN->unnu; i++)
		{
			while (1)
			{
				if (unpa[i].R == 0 || unpa[i].thrd == 0)
					break;
				else
					::Sleep(2);
			}
		}
		// 如果中途终止
		if (*(TAN->exit) == 1)
			break;
		// 构造隔离区
		R = (TIN_R*)B2->RemoveHead();
		while (R)
		{
			G = R->G;
			// 构造隔离区
			ce[0] = (G->bnd[0] + G->bnd[2]) / 2;
			ce[1] = (G->bnd[1] + G->bnd[3]) / 2;
			if (R->AXI == 0)	// 左右隔离
			{
				TIN_R *LR = new TIN_R;	// 左隔离区
				TIN_R *RR = new TIN_R;	// 右隔离区

				ce[0] = R->AXV;
				B1->AddTail(LR);	B1->AddTail(RR);
				ConRxn2(G, LR, RR, R->ASL, P, ce[0], ce[1]);
			}
			else				// 上下隔离
			{
				TIN_R *UR = new TIN_R;	// 上隔离区
				TIN_R *DR = new TIN_R;	// 下隔离区

				ce[1] = R->AXV;
				B1->AddTail(UR);	B1->AddTail(DR);
				ConRyn2(G, UR, DR, R->ASL, P, ce[0], ce[1]);
			}
			TAN->G->vpn += G->vpn;
			// 搜集三角形
			TIN_T *T = (TIN_T*)R->ATL->RemoveHead();
			while (T)
			{
				ATL->AddTail(T);
				T = (TIN_T*)R->ATL->RemoveHead();
			}
			// cid逆变换参数
			TCS_S32 sc = TCS_S32((G->bnd[0] - TAN->G->bnd[0]) / G->spn[0]);
			TCS_S32 sr = TCS_S32((G->bnd[1] - TAN->G->bnd[1]) / G->spn[1]);
			// 搜集静止边
			TIN_S *S = (TIN_S*)R->BSL->RemoveHead();
			while (S)
			{
				// cid逆变换
				S->sid[1] = InvCid2(S->sid[1], TAN->G->num[0], G->num[0], sc, sr);
				S->eid[1] = InvCid2(S->eid[1], TAN->G->num[0], G->num[0], sc, sr);
				BSL->PutData(S, S);
				S = (TIN_S*)R->BSL->RemoveNext();
			}
			// 释放隔离区内存资源
			delete R->ATL;	delete R->ASL;	delete R->BSL;
			delete []G->idx;	delete G;	delete R;
			// 如果中途终止
			if (*(TAN->exit) == 1)
				break;

			R = (TIN_R*)B2->RemoveHead();
		}
		// 如果中途终止
		if (*(TAN->exit) == 1)
			break;
	}
	// 处理中途终止
	TIN_R *R = (TIN_R*)B2->RemoveHead();
	while (R)
	{
		B1->AddTail(R);
		R = (TIN_R*)B2->RemoveHead();
	}
	// 释放内存资源
	delete []nven;	delete B2;
}

TCS_U32 WINAPI TanThread(LPVOID lpParam)
{
	TAN_PAR *TAN = (TAN_PAR*)lpParam;

	// 消除重合点
	RmvDup2(TAN);

	::_ftime_s(TAN->btm);

	*TAN->etm = *TAN->btm;
	// 生成构网分区
	if (*TAN->vpn >= 3)
		CnsTen2(TAN);

	::_ftime_s(TAN->etm);
	// 结束工作线程
	*(TAN->exit) = 1;
	TUN_PAR *unpa = (TUN_PAR*)TAN->unpa;
	for (long i = 0; i < TAN->unnu; i++)
	{
		::SetEvent(unpa[i].nven[0]);
		::WaitForSingleObject(unpa[i].thrd, INFINITE);
	}
	// 本线程结束
	TAN->thrd = 0;

	return 0;
}

/////////////////////////////////////////////////////////////////////////////
// CTcsTan2

CTcsTan2::CTcsTan2(TCS_V00)
{
	SYSTEM_INFO info;	::GetSystemInfo(&info);		m_cpu = info.dwNumberOfProcessors;

	m_P = 0;	m_nzn = 20000;		m_vpn = 0;

	m_G.idx = 0;	m_G.pid = 0;

	m_ATL = new CTcsBL;
	m_ARL = new CTcsBL;
	m_BSL = new CTcsHL(128, SIDE_EQLH, SIDE_HASH);

	// 分配线程参数
	long unnu = m_cpu;	if (unnu <= 2)	unnu = 2;
	TUN_PAR *unpa = new TUN_PAR[unnu];
	for (long i = 0; i < unnu; i++)
	{
		unpa[i].nven[0] = ::CreateEvent(0, 0, 0, 0);
		unpa[i].nven[1] = ::CreateEvent(0, 0, 1, 0);
		unpa[i].exit = &m_exit;
		unpa[i].thrd = 0;
		unpa[i].Z = 0;
		unpa[i].R = 0;
	}
	m_unpa = unpa;

	TAN_PAR *anpa = new TAN_PAR;
	anpa->exit = &m_exit;
	anpa->unnu = unnu;
	anpa->unpa = unpa;
	anpa->thrd = 0;

	anpa->nzn = &m_nzn;
	anpa->vpn = &m_vpn;
	anpa->btm = &m_btm;
	anpa->etm = &m_etm;
	anpa->ATL = m_ATL;
	anpa->ARL = m_ARL;
	anpa->BSL = m_BSL;
	anpa->G = &m_G;
	m_anpa = anpa;
}

CTcsTan2::~CTcsTan2(TCS_V00)
{
	m_exit = 1;
	TAN_PAR *anpa = (TAN_PAR*)m_anpa;
	// 等待线程结束
	::WaitForSingleObject(anpa->thrd, INFINITE);
	// 释放资源
	TUN_PAR *unpa = (TUN_PAR*)m_unpa;
	for (long i = 0; i < anpa->unnu; i++)
	{
		::CloseHandle(unpa[i].nven[0]);
		::CloseHandle(unpa[i].nven[1]);
	}
	delete []unpa;	delete anpa;

	Free();

	delete m_ATL;	delete m_ARL;	delete m_BSL;
}

TCS_V00 CTcsTan2::Free(TCS_V00)
{
	if (m_G.idx)
	{
		long n = m_G.num[0] * m_G.num[1];
		for (long i = 0; i < n; i++)
		{
			if (m_G.idx[i].num > 0)
				delete []m_G.idx[i].pid;
		}
		delete []m_G.idx;
	}
	m_G.idx = 0;
	if (m_G.pid)
		delete []m_G.pid;
	m_G.pid = 0;
	// 释放链表内存
	TIN_T *T = (TIN_T*)m_ATL->RemoveHead();
	while (T)
	{
		delete T;
		T = (TIN_T*)m_ATL->RemoveHead();
	}
	TIN_S *S = (TIN_S*)m_BSL->RemoveHead();
	while (S)
	{
		delete S;
		S = (TIN_S*)m_BSL->RemoveNext();
	}
	TIN_R *R = (TIN_R*)m_ARL->RemoveHead();
	while (R)
	{
		CTcsBL *ATL = R->ATL;
		T = (TIN_T*)ATL->RemoveHead();
		while (T)
		{
			delete T;
			T = (TIN_T*)ATL->RemoveHead();
		}
		delete ATL;
		CTcsHL *ASL = R->ASL;
		S = (TIN_S*)ASL->RemoveHead();
		while (S)
		{
			delete S;
			S = (TIN_S*)ASL->RemoveNext();
		}
		delete ASL;
		CTcsHL *BSL = R->BSL;
		S = (TIN_S*)BSL->RemoveHead();
		while (S)
		{
			delete S;
			S = (TIN_S*)BSL->RemoveNext();
		}
		delete BSL;
		delete R;

		R = (TIN_R*)m_ARL->RemoveHead();
	}
}

TCS_S32 CTcsTan2::Make(TCS_F64 *P, TCS_S32 N, TCS_F64 dup)
{
	if (P == 0 || N <= 0 || dup < 0)
		return 0;
	// 判断是否存在运行任务
	TAN_PAR *anpa = (TAN_PAR*)m_anpa;
	if (anpa->thrd != 0)
		return 0;

	// 释放资源
	Free();
	// 初始化变量
	m_vpn = N;
	m_G.vpn = 0;
	m_G.pid = new TCS_S32[N];
	m_P = P;	m_N = N;	m_G.dup = dup;
	::memset(m_G.pid, 0, sizeof(TCS_S32) * N);
	ComBnd2(&m_G, m_P, m_N);	CnsIdx2(&m_G, m_P, m_N, 4);
	// 创建线程
	m_exit = 0;
	TUN_PAR *unpa = (TUN_PAR*)m_unpa;
	// 创建工作线程
	for (long i = 0; i < anpa->unnu; i++)
	{
		::SetEvent(unpa[i].nven[1]);
		::ResetEvent(unpa[i].nven[0]);
		unpa[i].Z = 0;	unpa[i].R = 0;
		unpa[i].thrd = ::CreateThread(0, 0, TunThread, &unpa[i], 0, 0);
	}
	// 创建控制线程
	anpa->P = P;
	anpa->thrd = ::CreateThread(0, 0, TanThread, anpa, 0, 0);

	return 1;
}

TCS_F64 CTcsTan2::Time(TCS_V00)
{
	TCS_TMB etm;	::_ftime_s(&etm);

	TCS_F64 dtm = ::difftime(m_etm.time, m_btm.time) + (m_etm.millitm - m_btm.millitm);
	if (dtm == 0)
		return (::difftime(etm.time, m_btm.time) + (etm.millitm - m_btm.millitm) / 1000.0);
	else
		return (::difftime(m_etm.time, m_btm.time) + (m_etm.millitm - m_btm.millitm) / 1000.0);
}

TCS_V00 CTcsTan2::Stop(TCS_V00)
{
	m_exit = 1;
	TAN_PAR *anpa = (TAN_PAR*)m_anpa;
	// 等待线程结束
	::WaitForSingleObject(anpa->thrd, INFINITE);
}

TCS_S32 CTcsTan2::Done(TCS_V00)
{
	TAN_PAR *anpa = (TAN_PAR*)m_anpa;
	if (anpa->thrd == 0)
		return 1;

	return 0;
}

TCS_V00 CTcsTan2::Proc(TCS_S32 &PN, TCS_S32 &TN)
{
	PN = TN = 0;
	// 统计完成点数和三角形数
	PN = m_G.vpn;
	TN = m_ATL->GetSize();

	TIN_R *R = (TIN_R*)m_ARL->GetHead();
	while (R)
	{
		PN += R->G->vpn;
		TN += R->ATL->GetSize();
		R = (TIN_R*)m_ARL->GetNext();
	}
}

TIN_G *CTcsTan2::TinG(TCS_V00)
{
	return &m_G;
}

TIN_T *CTcsTan2::TinT(TCS_S32 &TN)
{
	// 统计三角形数目
	TN = m_ATL->GetSize();
	TIN_R *R = (TIN_R*)m_ARL->GetHead();
	while (R)
	{
		TN += R->ATL->GetSize();
		R = (TIN_R*)m_ARL->GetNext();
	}
	// 获得三角形数据
	TIN_T *T = new TIN_T[TN];

	TCS_S32 i = 0;
	TIN_T *E = (TIN_T*)m_ATL->GetHead();
	while (E)
	{
		T[i] = *E;	i++;
		E = (TIN_T*)m_ATL->GetNext();
	}
	R = (TIN_R*)m_ARL->GetHead();
	while (R)
	{
		E = (TIN_T*)R->ATL->GetHead();
		while (E)
		{
			T[i] = *E;	i++;
			E = (TIN_T*)R->ATL->GetNext();
		}
		R = (TIN_R*)m_ARL->GetNext();
	}

	return T;
}

TIN_R **CTcsTan2::TinR(TCS_S32 &RN)
{
	RN = m_ARL->GetSize();
	if (RN == 0)
		return 0;
	TIN_R **T = new TIN_R*[RN];

	TCS_S32 i = 0;
	TIN_R *R = (TIN_R*)m_ARL->GetHead();
	while (R)
	{
		T[i] = R;	i++;
		R = (TIN_R*)m_ARL->GetNext();
	}

	return T;
}