某高斯消元和某2-sat

poj3678
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
#include<string>
#include<iostream>
#include<stack>
#include<set>
#include<vector>
#include<map>
#include<queue>
using namespace std;
int dfn[2005],low[2005],mstack[2005],top;
bool instack[2005];
int head[2005],next[4000005],to[4000005],tot;
int dfsnum,cnt,color[2005],n,m,shu1,shu2,shu3;
char cha[10];
void tarjan(int x)
{
	dfn[x]=low[x]=++dfsnum;
	mstack[++top]=x;
	instack[x]=true;
	for(int i=head[x];i;i=next[i])
	{
		int y=to[i];
		if(!dfn[y])
		{
			tarjan(y);
			low[x]=min(low[x],low[y]);
		}
		else if(instack[y]) low[x]=min(low[x],dfn[y]);
	}
	if(dfn[x]==low[x])
	{
		cnt++;
		while(1)
		{
			int shu=mstack[top];
			top--;
			color[shu]=cnt;
			instack[shu]=false;
			if(shu==x)
			break;
		}
	}
}
void add(int x,int y)
{
	to[++tot]=y;
	next[tot]=head[x];
	head[x]=tot;
}
int main()
{
	cin>>n>>m;
	for(int i=1;i<=m;i++)
	{
		scanf("%d %d %d %s",&shu1,&shu2,&shu3,cha);
		shu1+=1;
		shu2+=1;
		if(cha[0]=='A')
		{
			if(shu3)
			{
				add(shu1*2,shu2*2);
				add(shu2*2,shu1*2);
				add(shu1*2-1,shu1*2);
				add(shu2*2-1,shu2*2);
			}
			else
			{
				add(shu1*2,shu2*2-1);
				add(shu2*2,shu1*2-1);
			}
		}
		else if(cha[0]=='O')
		{
			if(shu3)
			{
				add(shu2*2-1,shu1*2);
				add(shu1*2-1,shu2*2);
			}
			else
			{
				add(shu1*2-1,shu2*2-1);
				add(shu2*2-1,shu1*2-1);
				add(shu1*2,shu1*2-1);
				add(shu2*2,shu2*2-1);
			}
		}
		else if(cha[0]=='X')
		{
			if(shu3)
			{
				add(shu1*2,shu2*2-1);
				add(shu2*2,shu1*2-1);
				add(shu1*2-1,shu2*2);
				add(shu2*2-1,shu1*2);
			}
			else 
			{
				add(shu1*2,shu2*2);
				add(shu2*2,shu1*2);
				add(shu1*2-1,shu2*2-1);
				add(shu2*2-1,shu1*2-1);
			}
		}
	}
	for(int i=1;i<=2*n;i++)
	if(!dfn[i]) 
	tarjan(i);
	for(int i=1;i<=2*n;i+=2)
	{
		if(color[i]==color[i+1])
		{
			printf("NO");
			return 0;
		}
	}
	printf("YES");
	return 0;
}
/*
2 2
0 1 1 AND
0 1 1 XOR
*/
poj3207
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
#include<string>
#include<iostream>
#include<stack>
#include<set>
#include<vector>
#include<map>
#include<queue>
using namespace std;
int dfn[2005],low[2005],mstack[2005],top;
bool instack[2005];
int head[2005],next[4000005],to[4000005],tot,a[1005],b[1005];
int dfsnum,cnt,color[2005],n,m,shu1,shu2,shu3;
void add(int x,int y)
{
	to[++tot]=y;
	next[tot]=head[x];
	head[x]=tot;
}
void tarjan(int x)
{
	dfn[x]=low[x]=++dfsnum;
	mstack[++top]=x;
	instack[x]=true;
	for(int i=head[x];i;i=next[i])
	{
		int y=to[i];
		if(!dfn[y])
		{
			tarjan(y);
			low[x]=min(low[x],low[y]);
		}
		else if(instack[y]) low[x]=min(low[x],dfn[y]);
	}
	if(dfn[x]==low[x])
	{
		cnt++;
		while(1)
		{
			int shu=mstack[top];
			top--;
			color[shu]=cnt;
			instack[shu]=false;
			if(shu==x)
			break;
		}
	}
}
int main()
{
	cin>>n>>m;
	for(int i=1;i<=m;i++)
	{
		scanf("%d%d",&a[i],&b[i]);
		if(a[i]>b[i]) 
		swap(a[i],b[i]);
	}
	for(int i=1;i<=m;i++)
	for(int j=i+1;j<=m;j++)
	{
		if((a[i]<a[j]&&a[j]<b[i]&&b[i]<b[j])||(a[j]<a[i]&&a[i]<b[j]&&b[j]<b[i]))
		{
			add(2*i,2*j-1);
			add(2*j,2*i-1);
			add(2*j-1,2*i);
			add(2*i-1,2*j);
		}
	}
	for(int i=1;i<=2*m;i++)
	if(!dfn[i])
	tarjan(i);
	for(int i=1;i<=m*2;i+=2)
	if(color[i]==color[i+1])
	{
		printf("the evil panda is lying again");
		return 0;
	}
	printf("panda is telling the truth...");
	return 0;
}

以及poj1201

需要注意的是要建一条反向边表示sum[i]-sum[i-1]<=1

也就是add(i+1,i,-1)

进行限制

昨天的圆桌骑士和那道题终于找到了错误原因

当前弧优化应该hh[x]=i而不是next[i]

当前这条边可能流量并没有流完

今天又欠了一道题

WA掉的poj2723

明天继续写

寒假作业怕是写不完了

凉凉

#include<iostream> #include <fstream> #include <vector> #include <cmath> #include<sstream> #include <iomanip> #include <string> using namespace std; //矩阵转置 template<typename T, int x, int y> void transposition(const T(&a)[x][y], T(&at)[y][x]) { for (int i = 0; i < x; i++) { for (int j = 0; j < y; j++) { at[j][i] = a[i][j]; } } } //矩阵相乘 template<typename T, int x, int y, int z> void multiply(const T(&a)[x][y], const T(&b)[y][z], T(&c)[x][z]) { for (int i = 0; i < x; i++) { for (int j = 0; j < z; j++) { c[i][j] = 0; for (int k = 0; k < y; k++) { c[i][j] += a[i][k] * b[k][j]; } } } } //矩阵求逆 // 高斯-约旦元法求矩阵的逆 template<typename T, int size> bool inverse(const T(&a)[size][size], T(&inv)[size][size]) { // 创建增广矩阵 [A | I] T matrix[size][size * 2]; for (int i = 0; i < size; i++) { for (int j = 0; j < size; j++) { matrix[i][j] = a[i][j]; // 初始化单位矩阵部分 matrix[i][j + size] = (i == j) ? 1 : 0; } } // 高斯-约旦元 for (int p = 0; p < size; p++) { // 寻找主元 int max_row = p; for (int row = p + 1; row < size; row++) { if (abs(matrix[row][p]) > abs(matrix[max_row][p])) { max_row = row; } } // 交换 if (max_row != p) { for (int col = 0; col < 2 * size; col++) { swap(matrix[p][col], matrix[max_row][col]); } } // 检查是否为奇异矩阵 if (matrix[p][p] == 0) { cout << "矩阵是奇异的,无法求逆。" << endl; return false; } // 归一化主元 T pivot = matrix[p][p]; for (int col = 0; col < 2 * size; col++) { matrix[p][col] /= pivot; } // 元其他 for (int row = 0; row < size; row++) { if (row != p) { T factor = matrix[row][p]; for (int col = 0; col < 2 * size; col++) { matrix[row][col] -= factor * matrix[p][col]; } } } } // 提取逆矩阵部分 for (int i = 0; i < size; i++) { for (int j = 0; j < size; j++) { inv[i][j] = matrix[i][j + size]; } } return true; } // 常量定义(无法改变) const double c = 299792458.0;//光速 const double threshold = 0.0001;//迭代阈值 const int maxiteration = 100;//最大迭代次数 // 卫星数据结构体 struct SatelliteData { string ID;//卫星ID double x, y, z;// 卫星坐标 double range;// 伪距观测值 double variance;// 中误差 }; // 定义历元数据结构体,存储一个历元内的所有卫星数据 struct EpochData { vector<SatelliteData> satellites; //该历元内的卫星数据数组 }; int main() { // 打开数据文件 ifstream file("D:\\桌面\\数据\\卫星导航\\实验二\\expt2024.data"); if (!file.is_open()) { cout << "无法打开数据文件expt2024.data" << endl; return 1; } vector<EpochData> allEpochs; // 存储所有历元的数据 EpochData currentEpoch; // 当前处理的历元数据 string line; int satelliteCount = 0; // 当前历元已读取的卫星数 int epochCount = 0;//历元数 // 逐读取文件内容 while (getline(file, line) && epochCount < 5)//只计算了5个历元 { // 处理历元起始(以">"开头的) if (!line.empty() && line[0] == '>') { // 保存上一历元的数据 if (satelliteCount > 0) { allEpochs.push_back(currentEpoch); currentEpoch.satellites.clear(); satelliteCount = 0; epochCount++; } continue; // 跳过历元头,不读取具体数据 } // 处理卫星数据 else if (line.find("G") == 0) { // 卫星ID以"G"开头 istringstream iss(line); SatelliteData satData; // 按顺序读取卫星数据的各个字段 iss >> satData.ID >> satData.x >> satData.y >> satData.z >> satData.range >> satData.variance; // 将当前卫星数据添加到当前历元 currentEpoch.satellites.push_back(satData); satelliteCount++; } } file.close(); // 输出验证读取结果 cout << "文件读取成功,成功读取 " << allEpochs.size() << " 个历元的数据" << endl; //对5个历元进遍历计算接收机坐标 for (size_t i = 0; i < allEpochs.size(); i++) { vector<double> x, y, z; vector<double> p; vector<double> o; vector<double> P; //cout << "\n=== 历元 " << i + 1 << " ===" << endl; //cout << "包含 " << allEpochs[i].satellites.size() << " 颗卫星" << endl; for (size_t j = 0; j < allEpochs[i].satellites.size(); j++) { //提取数据为数组 const SatelliteData& sat = allEpochs[i].satellites[j]; x.push_back(sat.x); y.push_back(sat.y); z.push_back(sat.z); p.push_back(sat.range); o.push_back(sat.variance); } double X0[4][1] = { 0 }; // 初始坐标(X,Y,Z)和接收机钟差 double dX[4][1] = { 0 }; // 坐标增量和钟差增量 double X[4][1] = { 0 }; // 最终坐标和钟差 int iteration = 0; // 迭代次数 //迭代计算 while (iteration < maxiteration) { iteration++; double B[9][4] = { 0 }; // 系数矩阵B,每个历元9颗卫星 double L[9][1] = { 0 }; // 常数项L double P[9][9] = { 0 }; // 权矩阵P double BT[4][9] = { 0 };// 矩阵B的转置 double temp[4][9] = { 0 };//临时储存BT*P double N[4][4] = { 0 };//定义N=BT*P*B double invN[4][4] = { 0 };//定义矩阵N的逆矩阵invN double W[4][1] = { 0 };//定义矩阵W=BT*P*L int satelliteNum = 9;//每个历元中计算的卫星数 // 构建系数矩阵B、常数项L和权矩阵P for (int j = 0; j < satelliteNum; j++) { const SatelliteData& sat = allEpochs[i].satellites[j]; //读取每颗卫星 double x_s = sat.x; double y_s = sat.y; double z_s = sat.z; double range = sat.range; // 计算卫星到接收机的近似距离 double r0 = sqrt(pow(x_s - X0[0][0], 2) + pow(y_s - X0[1][0], 2) + pow(z_s - X0[2][0], 2)); // 计算系数矩阵B的元素 B[j][0] = -(x_s - X0[0][0]) / r0; B[j][1] = -(y_s - X0[1][0]) / r0; B[j][2] = -(z_s - X0[2][0]) / r0; B[j][3] = 1; // 接收机钟差系数 // 构建常数项L(伪距观测值已进卫星钟差改正,忽略电离层和对流层延迟) L[j][1] = range - r0; //计算权阵P double sigma = sat.variance; P[j][j] = 1.0 / (sigma * sigma); } //转置矩阵B transposition(B, BT); //计算BT*P multiply(BT, P, temp); multiply(temp, B, N); //对矩阵N求逆 inverse(N, invN); //计算BT*P*L multiply(temp, L, W); //计算dX multiply(N, W, dX); //更新接收机坐标值和接收机钟差 for (int k = 0; k < 4; k++) { X0[k][0] += dX[k][0]; } if (dX[0][0] < threshold && dX[1][0] < threshold && dX[2][0] < threshold) { cout << "计算第" << i + 1 << "个历元的数据一共迭代了" << iteration << "次,计算出的接收机坐标和接收机钟差三线表如下:" << endl; cout << "------------------------------------------------" << endl; cout << "----X----|----Y----|----Z----|----接收机钟差----" << endl; cout << "------------------------------------------------" << endl; cout << X0[0][0] << " | " << X0[1][0] << " | " << X0[2][0] << " | " << X0[3][0] << endl; cout << "------------------------------------------------" << endl; break; } } } return 0; }摄影测量单点定位程序设计
06-05
评论
成就一亿技术人!
拼手气红包6.0元
还能输入1000个字符
 
红包 添加红包
表情包 插入表情
 条评论被折叠 查看
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值