USACO Feed Ratios

1、这是98年的一道ACM的Final题目，可以暴力枚举，但是我想试试用高斯消元求非负整数解，所以就写了一个比较长的代码（当然运行时间很短）。

2、要注意，由于是求非负整数解，在两行相减的时候，其第一个数不一定是成倍数关系，所以要利用最小公倍数。注意此时循环应从0到n，否则会漏掉一些地方。计算后再将整行除以最大公因数，防止数字过大导致溢出。

3、最后还是要吐槽下：谁说USACO Training是给初学者用的题库啊啊啊，一会儿IOI一会儿Final我真受不了啊。。。

/*
ID:mrxy564
PROG:ratios
LANG:C++
*/
#include<cstdio>
#include<cmath>
#include<algorithm>
using namespace std;
int G[3][4],b[3],ans[3],num;
bool l[3];
int gcd(int a,int b){
    return b==0?a:gcd(b,a%b);
}
int lcm(int a,int b){
    return a*b/gcd(a,b);
}
int Gauss(int G[][4],bool l[],int ans[],const int &n){
    int res=0,r=0,a[3][4];
 for(int i=0;i<n;i++)
  for(int j=0;j<=n;j++)
   a[i][j]=G[i][j];
 for(int i=0;i<n;i++)
  l[i]=false;
 for(int i=0;i<n;i++){
  for(int j=r;j<n;j++)
   if(a[j][i]!=0){
       for(int k=i;k<=n;k++)
     swap(a[j][k],a[r][k]);
    break;
   }
 if(a[r][i]==0){
    res++;
    continue;
 }
 for(int j=0;j<n;j++)
  if(j!=r&&a[j][i]!=0){
   int temp=lcm(a[j][i],a[r][i]);
   int x1=temp/a[j][i],x2=temp/a[r][i];
   for(int k=0;k<=n;k++)
    a[j][k]=a[j][k]*x1-a[r][k]*x2;
   int g=a[j][i];
   for(int k=0;k<=n;k++)
    g=gcd(g,a[j][k]);
   for(int k=0;k<=n;k++)
    a[j][k]/=g;
  }
     l[i]=true;r++;
 }
 for(int i=0;i<n;i++)
  if(l[i])
   for(int j=0;j<n;j++)
    if(a[j][i]!=0){
     if(a[j][n]%a[j][i]==0)
        ans[i]=a[j][n]/a[j][i];
     else return -1;
                                        if(ans[i]<0) return -1;
    }
 return res;
}
int main(){
 freopen("ratios.in","r",stdin);
 freopen("ratios.out","w",stdout);
 for(int j=0;j<3;j++){
  scanf("%d",&b[j]);
  G[j][3]=0;
 }
    for(int i=0;i<3;i++)
  for(int j=0;j<3;j++)
   scanf("%d",&G[j][i]);
        bool flag1=false;
 for(int x=1;x<100;x++){
     for(int j=0;j<3;j++)
   G[j][3]=b[j]*x;
  num=Gauss(G,l,ans,3);
  bool flag=true;
  if(num<0) flag=false;
  else{
                           for(int k=3;k>3-num;k--)
      if(G[k][3]) {flag=false;break;}
  }
  if(flag){
      printf("%d %d %d %d\n",ans[0],ans[1],ans[2],x);
                    flag1=true;
   break;
  }
 }
        if(!flag1) printf("NONE\n");
 return 0;
}
官方题解：

As there are only 101³ = 1,030,301 ways to do this, try them all and pick the best.

#include <stdio.h>

/* the goal ratio */
int goal[3];

/* the ratio of the feeds */
int ratios[3][3];

/* the best solution found so far */
int min;
int mloc[4]; /* amounts of ratio 1, 2, and 3, and the amount of ratio 4 prod */

/* the amount of each type of component in the feed */
int sum[3];

int main(int argc, char **argv)
 {
  FILE *fout, *fin;
  int lv, lv2, lv3; /* loop variables */
  int t, s; /* temporary variables */
  int gsum; /* the sum of the amounts of components in the goal mixture */

  if ((fin = fopen("ratios.in", "r")) == NULL)
   {
    perror ("fopen fin");
    exit(1);
   }
  if ((fout = fopen("ratios.out", "w")) == NULL)
   {
    perror ("fopen fout");
    exit(1);
   }

  /* read in data */
  fscanf (fin, "%d %d %d", &goal[0], &goal[1], &goal[2]);
  for (lv = 0; lv < 3; lv++)
    fscanf (fin, "%d %d %d", ratios[lv]+0, ratios[lv]+1, ratios[lv]+2);

  gsum = goal[0] + goal[1] + goal[2];

  min = 301; /* worst than possible = infinity */

  /* boundary case (this ensures gsum != 0) */
  if (gsum == 0)
   {
    fprintf (fout, "0 0 0 0\n");
    return 0;
   }

  for (lv = 0; lv < 100; lv++)
    for (lv2 = 0; lv2 < 100; lv2++)
     { /* for each amout of the first two types of mixtures */
      sum[0] = lv*ratios[0][0] + lv2*ratios[1][0];
      sum[1] = lv*ratios[0][1] + lv2*ratios[1][1];
      sum[2] = lv*ratios[0][2] + lv2*ratios[1][2];

      if (lv + lv2 > min) break;

      for (lv3 = 0; lv3 < 100; lv3++)
       {
        s = lv + lv2 + lv3;
	if (s >= min) break; /* worse than we've found already */

        /* calculate a guess at the multiples of the goal we've obtained */
	/* use gsum so we don't have to worry about divide by zero */
        t = (sum[0] + sum[1] + sum[2]) / gsum;
	if (t != 0 && sum[0] == t*goal[0] && 
	        sum[1] == t*goal[1] && sum[2] == t*goal[2])
	 { /* we've found a better solution! */
	  /* update min */
	  min = s;

	  /* store the solution */
	  mloc[0] = lv;
	  mloc[1] = lv2;
	  mloc[2] = lv3;
	  mloc[3] = t;
	 }

        /* add another 'bucket' of feed #2 */
        sum[0] += ratios[2][0];
        sum[1] += ratios[2][1];
        sum[2] += ratios[2][2];
       }
     }
  if (min == 301) fprintf (fout, "NONE\n"); /* no solution found */
  else fprintf (fout, "%d %d %d %d\n", mloc[0], mloc[1], mloc[2], mloc[3]);
  return 0;
 }

Vlad Novakovski's Solution

When you combine multiples of mixtures, you can look at it as a multiplication of a matrix by a vector. For example, 8*(1:2:3) + 1*(3:7:1) + 5*(2:1:2) = (21:28:35) = 7*(3:4:5) can be seen as

[ 1 3 2 ]   [ 8 ]            
[ 2 7 1 ] * [ 1 ] = 7 [3 4 5]
[ 3 1 2 ]   [ 5 ]           

The matrix and the goal ratio vector (3:4:5 in this case) are given; what we have to find is the multiple vector (8:1:5 in this case) and the proportionality costant (7 here). This is like solving a system of linear equations. We can write it as

AX = kB.

Now, if we use Cramer's Rule, and let D = determinant of A, then

X_1 = k D_1 / D
X_2 = k D_2 / D
X_3 = k D_3 / D,

where D_1 is the determinant of the matrix A with the 1st column is replaced by B, D_2 is the determinant of the matrix A with the 2nd column is replaced by B, D_3 is the determinant of the matrix A with the 3rd column is replaced by B. (see a Linear Algebra textbook on why this works.) ,P> We are looking for integral solutions. If D = 0, no solutions. Otherwise, let k = D, and then X_1 = D_1, etc. If these values (X_1,X_2,X_3, _and_ k) all have a greatest common factor above 1, divide them all by that factor, as we are looking for the smallest possible solutions.

Now, if some of the numbers is greater than 100, we have not found a feasible solution, so we output `NONE'. Otherwise, the triple (X_1,X_2,X_3) is output, as well as the value k.