Shen Huang's Blog

/* 自适应辛普森 UVA 1356 
 * 二分求解高度 求解抛物线的长度是不是合适
 *
 * */
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <queue>
#include <set>
#include <vector>
using namespace std;

const int INF = ~0u>>1;
typedef pair <int,int> P;
#define MID(x,y) ((x+y)>>1)
#define iabs(x) ((x)>0?(x):-(x))
#define REP(i,a,b) for(int i=(a);i<(b);i++)
#define FOR(i,a,b) for(int i=(a);i<=(b);i++)
#define pb push_back
#define mp make_pair
#define print() cout<<"--------"<<endl

// 这里为了方便，把a声明成全局的。
// 这不是一个好的编程习惯，但在本题中却可以提高代码的可读性
double a; 

//simpson 公式用的函数 (这个需要根据题目而改变函数)
double F(double x){
	return sqrt(1 + 4*a*a*x*x);
}

//三点simpson法，这里要求F是一个全局函数
double simpson(double a,double b){
	double c = a + (b-a)/2;
	return (F(a) + 4*F(c) + F(b)) * (b-a) / 6;
}

//自适应simpson公式(递归)。已知整个区间[a,b]上三点的simpson值A
double asr(double a,double b,double eps,double A){
	double c = a + (b-a)/2;
	double L = simpson(a,c),R = simpson(c,b);
	if (fabs(L+R-A) <= 15*eps) return L+R+(L+R-A)/15.0;
	return asr(a,c,eps/2,L) + asr(c,b,eps/2,R);
}

//自适应simpson公式(主过程)
double asr(double a,double b,double eps){
	return asr(a,b,eps,simpson(a,b));
}

//下面这个函数依题面而异，本题求解simpson求解宽度为w,高度为h的抛物线长度
double solve(double w,double h){
	a = 4.0 * h / (w*w); //修改全局变量a，从而改变F的行为
	return asr(0,w/2,1e-5)*2;
}

int main(){
	int t,cas = 0;
	scanf("%d", &t);
	while (cas++<t){
		int D,H,B,L;
		scanf("%d%d%d%d", &D, &H, &B, &L);
		int n = (B+D-1) / D;
		double D1 = (double) B / n;
		double L1 = (double) L / n;
		double x = 0, y = H;
		while(y-x > 1e-5){//二分求解高度
			double m = x + (y-x) / 2;
			if(solve(D1,m) < L1) x = m;
			else y = m;
		}
		if(cas > 1) puts("");
		printf("Case %d:n%.2fn", cas, H - x);
	}

	return 0;

}

/* BGstep HDU 2815
 * 求解a^x == c (mod p) 的 0 <= x < c
 * */
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <queue>
#include <set>
#include <vector>
using namespace std;

const int INF = ~0u>>1;
typedef pair <int,int> P;
typedef long long llg;
#define MID(x,y) ((x+y)>>1)
#define iabs(x) ((x)>0?(x):-(x))
#define REP(i,a,b) for(int i=(a);i<(b);i++)
#define FOR(i,a,b) for(int i=(a);i<=(b);i++)
#define pb push_back
#define mp make_pair
#define print() cout<<"--------"<<endl
const int N = 100003;
struct Baby_step_Giant_step{
	struct edge{int key,times;edge*next;}g[N],*ls[N];
	int e;

	void add(int x,int times){
		int hash = x % N;
		for (edge*t=ls[hash];t;t=t->next)
			if(t->key == x) return ;
		g[++e].key = x; 
		g[e].times = times;
		g[e].next = ls[hash];
		ls[hash] = &g[e];
	}

	int In(int x){
		int hash = x % N;
		for (edge*t=ls[hash];t;t=t->next)
			if(t->key == x) return t->times;
		return -1;
	}

	int Power_mod(int a,int b,int mod){
		llg t = a, res = 1;
		while (b){
			if (b&1) res = res * t % mod;
			t = t * t % mod;
			b >>= 1;
		}
		return (int)(res);
	}

	llg exgcd(int a,int b,llg &x,llg &y){
		if (b == 0){
			x = 1; y = 0; return a; 
		}
		llg g = exgcd(b,a%b,x,y);
		llg tx = x,ty = y;
		x = ty; y = tx - a/b*ty;
		return g;
	}

	int solve(int a,int c,int p){
		if (c >= p) return -1;
		int i,j;
		int x = 1; int t = 1;
		int s = (int)ceil(sqrt((double)p));
		int res,cnt = 0;
		llg tx,ty,g;
		for (int i = 0; i <= 50; i++)
			if (Power_mod(a,i,p) == c) 
				return i;

		while ((g=exgcd(a,p,tx,ty)) != 1){
			if (c%g) return -1;
			c /= g;
			p /= g;
			x = (llg)(a)/g*x%p;
			cnt ++;
		} 

		e = 0;
		memset(ls,0,sizeof(ls));
		for (i = 0; i < s; i++) add(Power_mod(a,i,p),i);

		for (int step = Power_mod(a,s,p),i = cnt; i < p;
				i += s,x = (llg)(x) * step % p){
			g = exgcd(x,p,tx,ty);
			tx = tx * c % p;
			if (tx < 0) tx += p;
			res = In((int)tx);
			if (res == -1) continue;
			return res + i;
		}
		return -1;
	}
}Number_Theory;

int main(){
	int a,p,c;
	while(~scanf("%d%d%d", &a, &p, &c)){
		int ans = Number_Theory.solve(a,c,p);
		if (ans == -1) puts("Orz,I can&rsquo;t find D!");
		else printf("%dn", ans);
	}
	return 0;
}

#include <cstdio>
#include <cmath>
/**
Simplex C(n+m)(n)
maximize:
    c[1]*x[1]+c[2]*x[2]+...+c[n]*x[n]=ans
subject to
    a[1,1]*x[1]+a[1,2]*x[2]+...a[1,n]*x[n] <= rhs[1]
    a[2,1]*x[1]+a[2,2]*x[2]+...a[2,n]*x[n] <= rhs[2]
    ......
    a[m,1]*x[1]+a[m,2]*x[2]+...a[m,n]*x[n] <= rhs[m]
限制:
    传入的矩阵必须是标准形式的, 即目标函数要最大化；
	约束不等式均为<= ；xi为非负数(>=0).
simplex返回参数：
    OPTIMAL             有唯一最优解
    UNBOUNDED           最优值无边界
    FEASIBLE            有可行解
    INFEASIBLE          无解
n为元素个数,m为约束个数
线性规划：   
    max c[]*x;
    a[][]<=rhs[];
ans即为结果,x[]为一组解(最优解or可行解)
**/
const double eps = 1e-8;
const double inf = 0x3f3f3f3f;

#define OPTIMAL -1          //表示有唯一的最优基本可行解
#define UNBOUNDED -2        //表示目标函数的最大值无边界
#define FEASIBLE -3         //表示有可行解
#define INFEASIBLE -4       //表示无解
#define PIVOT_OK 1          //还可以松弛
#define maxn 1000

struct LinearProgramming{
    int basic[maxn], row[maxn], col[maxn];
    double c0[maxn];

    double dcmp(double x){
        if (x > eps)    return 1;
        else if (x < -eps)   return -1;
        return 0;
    }
    void init(int n, int m, double c[], double a[maxn][maxn], double rhs[], double &ans) { //初始化
        for(int i = 0; i <= n+m; i++) {
            for(int j = 0; j <= n+m; j++)   a[i][j]=0;
            basic[i]=0; row[i]=0; col[i]=0;
            c[i]=0; rhs[i]=0;
        }
        ans=0;
    }
    //转轴操作
    int Pivot(int n, int m, double c[], double a[maxn][maxn], double rhs[], int &i, int &j){
        double min = inf;
        int k = -1;
        for (j = 0; j <= n; j ++)
            if (!basic[j] && dcmp(c[j]) > 0)
                if (k < 0 || dcmp(c[j] - c[k]) > 0)     k = j;
        j = k;
        if (k < 0)  return OPTIMAL;
        for (k = -1, i = 1; i <= m; i ++)   if (dcmp(a[i][j]) > 0)
            if (dcmp(rhs[i] / a[i][j] - min) < 0){  min = rhs[i]/a[i][j]; k = i;  }
        i = k;
        if (k < 0)  return UNBOUNDED;
        else    return PIVOT_OK;
    }
    int PhaseII(int n, int m, double c[], double a[maxn][maxn], double rhs[], double &ans, int PivotIndex){
        int i, j, k, l; double tmp;
        while(k = Pivot(n, m, c, a, rhs, i, j), k == PIVOT_OK || PivotIndex){
            if (PivotIndex){    i = PivotIndex; j = PivotIndex = 0; }
            basic[row[i]] = 0;  col[row[i]] = 0;    basic[j] = 1;   col[j] = i;     row[i] = j;
            tmp = a[i][j];
            for (k = 0; k <= n; k ++)   a[i][k] /= tmp;
            rhs[i] /= tmp;
            for (k = 1; k <= m; k ++)
                if (k != i && dcmp(a[k][j])){
                    tmp = -a[k][j];
                    for (l = 0; l <= n; l ++)   a[k][l] +=  tmp*a[i][l];
                    rhs[k] += tmp*rhs[i];
                }
            tmp = -c[j];
            for (l = 0; l <= n; l ++)   c[l] += a[i][l]*tmp;
            ans -= tmp * rhs[i];
        }
        return k;
    }
    int PhaseI(int n, int m, double c[], double a[maxn][maxn], double rhs[], double &ans){
        int i, j, k = -1;
        double tmp, min = 0, ans0 = 0;
        for (i = 1; i <= m; i ++)
            if (dcmp(rhs[i]-min) < 0){min = rhs[i]; k = i;}
        if (k < 0)  return FEASIBLE;
        for (i = 1; i <= m; i ++)   a[i][0] = -1;
        for (j = 1; j <= n; j ++)   c0[j] = 0;
        c0[0] = -1;
        PhaseII(n, m, c0, a, rhs, ans0, k);
        if (dcmp(ans0) < 0) return INFEASIBLE;
        for (i = 1; i <= m; i ++)   a[i][0] = 0;
        for (j = 1; j <= n; j ++)
            if (dcmp(c[j]) && basic[j]){
                tmp = c[j];
                ans += rhs[col[j]] * tmp;
                for (i = 0; i <= n; i ++)   c[i] -= tmp*a[col[j]][i];
            }
        return FEASIBLE;
    }
    //standard form
    //n:原变量个数   m:原约束条件个数
    //c:目标函数系数向量-[1~n],c[0] = 0;
    //a:约束条件系数矩阵-[1~m][1~n]      rhs:约束条件不等式右边常数列向量-[1~m]
    //ans:最优值   x:最优解||可行解向量-[1~n]
    int simplex(int n, int m, double c[], double a[maxn][maxn], double rhs[], double &ans, double x[]){
        int i, j, k;
        //标准形式变松弛形式
        for (i = 1; i <= m; i ++){
            for (j = n+1; j <= n+m; j ++)   a[i][j] = 0;
            a[i][n+i] = 1;  a[i][0] = 0;
            row[i] = n+i;   col[n+i] = i;
        }
        k = PhaseI(n+m, m, c, a, rhs, ans);
        if (k == INFEASIBLE)    return k;
        k = PhaseII(n+m, m, c, a, rhs, ans, 0);
        for (j = 0; j <= n+m; j ++) x[j] = 0;
        for (i = 1; i <= m; i ++)   x[row[i]] = rhs[i];
        return k;
    }
}ps;    //Primal Simplex

int n,m;
double c[maxn], ans, a[maxn][maxn], rhs[maxn], x[maxn];

int main(){
    while(scanf("%d %d", &n, &m) != EOF){
        double ans;
        ps.init(n, m, c, a, rhs, ans);
        for (int i = 1; i <= n; i ++)
            scanf("%lf", &c[i]);
        for (int i = 1; i <= m; i ++){
            for (int j = 1; j <= n; j ++){
                scanf("%lf", &a[i][j]);
            }
            scanf("%lf", &rhs[i]);
        }
        int hehe=ps.simplex(n,m,c,a,rhs,ans,x);
        printf("Nasa can spend %.0f taka.n", ceil(m*ans));
    }
    return 0;
}

/* POJ 2184 DP 背包 
 * 在于有负数的情况，平移可以解决,平移100000，这样原来的0就是现在的100000
 * 在一般的01背包压缩空间的时候，体积的遍历是从大到小，
 * 因为dp[v]=max(dp[v],dp[v-c[i]]+w[i])，
 * 当前的dp[v]只取决于比自己小的dp[v-c[i]]，
 * 所以从大到小遍历时每次dp[v-c[i]]和dp[v]都是上一次的状态。
 * 如果体积为负v-c[i]>v，从大到小遍历dp[v-c[i]]是当前物品的状态，
 * 不是上一个，这样就会出错，解决的办法是从小到大遍历。
 *
 * 还有边界的部分不需要怎么考虑因为，数据比较小，
 * 按理不会出现达到-100000还能变到不是负数的情况所以你的dp[0]或者dp[p[i].s]
 * 开始没有多少区别
 * */
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <queue>
#include <set>
#include <vector>
#include <climits>
using namespace std;

const int INF = ~0u>>1;
typedef pair <int,int> P;
#define MID(x,y) ((x+y)>>1)
#define iabs(x) ((x)>0?(x):-(x))
#define REP(i,a,b) for(int i=(a);i<(b);i++)
#define FOR(i,a,b) for(int i=(a);i<=(b);i++)
#define pb push_back
#define mp make_pair
#define print() cout<<"--------"<<endl
const int maxn = 100005;
int dp[maxn * 2];
int n;
struct node{
	int s,f;
}p[200];
int main(){
	while (~scanf("%d", &n)){
		for (int i = 0; i < 2 * maxn; i++)
			dp[i] = INT_MIN;
		for (int i = 0; i < n; i++)
			scanf("%d%d", &p[i].s, &p[i].f);
		dp[100000] = 0;
		for (int i = 0; i < n; i++){
			if(p[i].s >= 0){
				for (int j = 2 * maxn - 1; j >= p[i].s; j--){
					if (dp[j-p[i].s] > INT_MIN)
						dp[j] = max(dp[j],dp[j-p[i].s]+p[i].f);
				}
			}else {
				for (int j = 0; j - p[i].s <= 200000; j++){
					if (dp[j-p[i].s] > INT_MIN)
						dp[j] = max(dp[j],dp[j-p[i].s]+p[i].f);
				}
			}
		}
		int ans = 0;
		for (int i = 100000; i < 2 * maxn; i++)
			if(dp[i] >= 0) // 注意f要>=0按照题意
				ans = max(ans,dp[i]+i-100000);
		printf("%dn", ans);
	}
	return 0;
}

/*
 *
 * */
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <queue>
#include <set>
#include <map>
#include <vector>
using namespace std;

const int INF = 10020000;
typedef pair <int,int> P;
#define MID(x,y) ((x+y)>>1)
#define iabs(x) ((x)>0?(x):-(x))
#define REP(i,a,b) for(int i=(a);i<(b);i++)
#define FOR(i,a,b) for(int i=(a);i<=(b);i++)
#define pb push_back
#define mp make_pair
#define print() cout<<"--------"<<endl
#define maxn 220
int n,m,num;
map<string,int>maps;
int edge[maxn][maxn];
int query(string s)
{
	map<string,int>::iterator it;
	it = maps.find(s);
	if (it == maps.end()){
		maps[s] = num++;
		return maps[s];
	}else return maps[s];
}

void floyd()
{
	for (int k = 1; k <= n; k++)
		for (int i = 1; i <= n; i++)
			for (int j = 1; j <= n; j++)
			{
				edge[i][j] = max(edge[i][j],min(edge[i][k],edge[k][j]));
			}
}

int main(){
	int co = 0,w;
	char s1[100],s2[100];
	while (scanf("%d%d", &n, &m) && (m + n)){
		getchar();
		num = 1;
		for (int i = 1; i <= n; ++i)
			for (int j = 1; j <= n; j++)
				edge[i][j] = 0;
		for (int i = 1; i <= m; ++i)
		{
			scanf("%s%s%d", s1, s2, &w);getchar();
			int x = query(s1);
			int y = query(s2);
			edge[x][y] = w;
			edge[y][x] = w;
		}
		floyd();
		scanf("%s%s",s1,s2);
		int x = query(s1);
		int y = query(s2);
		printf("Scenario #%dn%d tonsnn", ++co, edge[x][y]);
	}
	return 0;
}

//首先求末尾零的个数的求法，就是求5的因子的个数，因为2×5得到一个0.
//这个函数函数求1-N因子个数的求法
int cali(int n,int i){
 int ret = 0;
 while(n){
  n /= i;
  ret += n;
 }
 return ret;
}

void solve(int n){
 if(n<=1) return ;
 for (int i = n; i > 0; i /= 5){
  int p = i/10,q = i%10;
  co3 += p + (int)(q>=3);
  co5 += p + (int)(q>=5);
  co7 += p + (int)(q>=7);
  co9 += p + (int)(q>=9);
 }
 co2 += n/2;
 solve(n/2);
}

inline int pmod(int n,int k){
    int ret = 1;
    while (k){
        if (k&1) ret = ret * n % 10;
        n = n * n %10;
        k >>= 1;
    }
    return ret;
}