29 Temmuz 2018

Euler Projesi 246. Soru

Elipsin Teğetleri

Elipsin bir tanımı: M merkezli r yarıçaplı bir c çemberi ve d(G,M)<r olacak şekilde bir G noktası verilmek üzere, c ve G'den eş uzaklıktaki noktaların geometrik yeri.

Aşağıda tanıma ilişkin bir gösterim verilmektedir:
M(-2000,1500) ve G(8000,1500) noktaları veriliyor. Ayrıca M merkezli ve 15000 yarıçaplı c çemberi veriliyor. c ve G'den eş uzaklıktaki noktaların geometrik yeri e elipsi olsun. Elipsin dışındaki bir P noktasından elipse iki t1 ve t2 teğeti çiziliyor. Teğetlerin elipse değme noktaları R ve S olsun.
Kaç farklı P örgü noktası için RPS açısı 45 dereceden büyüktür?

(Project Euler)

Cevap: 810834388

Python:
import math

a = 5
b = 9
c = 281250000
k = 1/2

s = 0

y = 0
while True:
    yy = y * y
    xx = (a*c*(1+k)+b*(c-a*yy)*(1-k)+
          2*math.sqrt(a*c*k*(a*c+b*yy*(b-a)*(1-k))))/(a*b*(1-k))
    if xx < 0:
        break
    xx2 = (c-b*yy)/a
    if xx2 >= 0:
        i = 2 * (int(math.sqrt(xx)) - int(math.sqrt(xx2)))
    else:
        i = 2 * int(math.sqrt(xx)) + 1
    s += 2 * i if y else i
    y += 1

print(s)
C#:
static void Euler246()
    {
      double Ta = 7500;
      double Te = 5000;

      // a**2, b**2 e**2
      double a2 = Ta * Ta;
      double e2 = Te * Te;
      double b2 = a2 - e2;

      // get points with y=0 and with x=0 with exactly 45 degree angle.
      // these 2 values (a,b) forms an assumed ellipse. 
      // (i assumed these ellipse presents the border, but this is not the case)
      double aE2 = E246_getXP2(a2, b2);
      double bE2 = E246_getXP2(b2, a2);

      C.WriteLine("AE2=" + aE2);
      C.WriteLine("BE2=" + bE2);

      // calculate number of lattince points inside inner ellipse
      ulong lOI = E246_latticePointsInE(Ta, Math.Sqrt(b2));

      // now calculate for each y an x which is the border. (outside the angle < 45 degree)
      long lCnt = 0;
      for (long y = 1; y <= bE2; y++)
      {
        long x;
        // calculate the x value for this y value
        double xAssumed = (Math.Sqrt(bE2 * bE2 - y * y)) * aE2 / bE2;
        for (x = (long)xAssumed; E246_check_angle(a2, b2, x, y); x++) ;
        x--;

        //C.WriteLine("y=" + y + ",xR=" + x + "  ( xAssumed= " + xAssumed + ")");
        lCnt += x;
      }
      lCnt *= 4; // for all 4 qudarants
      lCnt += 2 * (long)aE2;  // for x axis, without origin
      lCnt += 2 * (long)bE2;  // for y axis, without origin
      lCnt++;  // origin

      // present the result.
      C.WriteLine("lCnt=" + lCnt);
      C.WriteLine(lCnt - (long)lOI);
    }

    static double E246_getXP2(double a, double b)
    {
      double T = Math.Tan(45.0*Math.PI/(2.0*180));
      double T2 = T*T;
      double p = -2 * a - b / T2;
      double q = (a * a) + (a * b) / T2;

      double mpd2 = -p / 2.0;
      double x1 = mpd2 + Math.Sqrt(mpd2 * mpd2 - q);
      double xp = Math.Sqrt(x1);
      return xp;
    }

    // calculate number of lattice points inside ellipse
    static ulong E246_latticePointsInE(double a, double b)
    {
      ulong lCnt = 0;
      for (ulong i = 1; i < b; i++)
      {
        double x=(Math.Sqrt(b * b - i * i)) * a / b;
        //C.WriteLine("---" + x);
        lCnt += (ulong)x;
      }
      //C.WriteLine(lCnt);
      lCnt *= 2;
      lCnt += (ulong)a;
      lCnt *= 2;
      lCnt += (ulong)b;
      lCnt += (ulong)b;
      lCnt++;

      return lCnt;
    }

    // for given x at ellipse (a**2,b**2) calculate y value
    static double E246_y(double x, double a2, double b2)
    {
      return Math.Sqrt((a2 - x * x) * b2 / a2);
    }

    // for given ellipse (a**2, b**2) and given point p(xp,yp) calculate
    // the angle between the tangents and check against 45 degrees.
    // return true if >= 45degrees
    static bool E246_check_angle(double a2, double b2, long xP, long yP)
    {
      double x1, x2, y1, y2;

      if (yP != 0)
      {
        double k = -b2 / a2 * xP / yP;
        double k2 = k * k;
        double d = b2 / yP;
        double v1 = a2 * k * d / (b2 + a2 * k2);
        double v2 = Math.Sqrt((a2 * k2 + b2 - d * d) * a2 * b2);
        v2 /= b2 + a2 * k2;
        x1 = -v1 + v2;
        x2 = -v1 - v2;

        y1 = ((1 - x1 * xP / a2) * b2) / yP;
        y2 = ((1 - x2 * xP / a2) * b2) / yP;
      }
      else
      {
        x1 = a2 / xP;
        x2 = x1;
        y1 = E246_y(x1,a2,b2);
        y2 = -y1;
      }

      double vx = x1 - xP;
      double vy = y1 - yP;
      double l1=Math.Sqrt(vx * vx + vy * vy);
      vx = x2 - xP;
      vy = y2 - yP;
      double l2 = Math.Sqrt(vx * vx + vy * vy);
      vx = x1 - x2;
      vy = y1 - y2;
      double l3 = Math.Sqrt(vx * vx + vy * vy);

      //C.WriteLine("P1:" + x1 + "," + y1);
      //C.WriteLine("P2:" + x2 + "," + y2);
      //C.WriteLine("P3:" + xP + "," + yP);

      double angle=Math.Acos((l1 * l1 + l2 * l2 - l3 * l3) / (2 * l1 * l2)) * 180 / Math.PI;
      //C.WriteLine("angle="+angle);

      if (angle > 45)
        return true;



      return false;
    }

Java:
import java.applet.Applet;
import java.awt.Color;
import java.awt.Font;
import java.awt.Frame;
import java.awt.Graphics;
import java.awt.Image;
import java.awt.event.WindowAdapter;
import java.awt.event.WindowEvent;
import java.util.ArrayList;
import java.util.List;

public class Euler246 extends Applet implements Runnable {
    class Triplet {
        public int x, y, size;
        public Triplet(int x, int y, int size)
            { this.x = x; this.y = y; this.size = size; }
        public void drawOn(Graphics g, Color c) {
            g.setColor(c);
            g.fillRect((x+16000) / 64, (y+21600) / 64, size / 64, size / 64);
        }
        public void proliferate(List list) {
            int s2 = size / 2;
            list.add(new Triplet(x,      y,      s2));
            list.add(new Triplet(x,      y + s2, s2));
            list.add(new Triplet(x + s2, y,      s2));
            list.add(new Triplet(x + s2, y + s2, s2));
        }
    }
        
    private Image img;
    private Thread kicker;
    
    public void paint(Graphics g) {
        if (img != null) {
            g.drawImage(img, 0, 0, this);
        }
    }
    
    public void update(Graphics g) {
        paint(g);
    }
    
    int judge(int x, int y) {
        if ((long)x * x * 5+ (long)y * y * 9 < 281250000) return 0;

        double a = 56250000 - x * x;
        double b = x * y;
        double c = 31250000 - y * y;

        int f = (a + c < 0) ? -1 : 1;
        
        return (f * (25*a*x*x - 90*b*x*y + 81*c*y*y) < 0 ||
                 c*c + 6*a*c - 4*b*b + a*a < 0) ? 1 : 2;
    }
    
    public void run() {
        img = createImage(getWidth(), getHeight());
        Graphics g = img.getGraphics();
        List list = new ArrayList();
        Color colors[] = new Color[] {
                new Color(0xFFFFCC), new Color(0xCCCCFF), Color.white }; 
        
        for (int ix = -20480; ix < 20480; ix += 4096)
            for (int iy = -20480; iy < 20480; iy += 4096)
                list.add(new Triplet(ix, iy, 4096));
        
        int count = 0;
        Font f = new Font("Serif", Font.PLAIN, 24);
        g.setFont(f);

        while (list.size() > 0 && kicker != null) {
            List next = new ArrayList();
            for (Triplet t : list) {
                int j = judge(t.x, t.y);
                if (t.size == 1 ||
                        (j == judge(t.x,             t.y + t.size - 1) &&
                         j == judge(t.x + t.size - 1, t.y            ) &&
                         j == judge(t.x + t.size - 1, t.y + t.size - 1))) {
                    t.drawOn(g, colors[j]);
                    if (j == 1) count += t.size * t.size;
                } else t.proliferate(next);
            }
            g.setColor(Color.white);
            g.fillRect(0, 640, 700, 36);
            g.setColor(Color.black);
            String msg = String.format("count=%d, rest blocks =%d", count, next.size()); 
            g.drawString(msg, 8, 660);
            repaint();
            try { Thread.sleep(500); } catch (Exception e) { }
            list = next;
        }
    }
        
    public void start() {
        kicker = new Thread(this);
        kicker.start();
    }
    
    public void stop() {
        kicker = null;
    }
    
    public static void main(String args[]) {
        Frame fr = new Frame();
        Euler246 inst = new Euler246();
        inst.setVisible(true);
        fr.add(inst);
        fr.setSize(500, 700);
        fr.setVisible(true);
        fr.addWindowListener(new WindowAdapter(){
            public void windowClosing(WindowEvent e) { System.exit(0);  }
        });
        inst.start();
    }
}