Página 3 de 4

MensagemEnviado: Sexta Jan 04, 2008 8:36 am
por AlexandreH
eu nao vejo programação, mas achei o código em Linguajem Java do procedimento que permite calcular a raiz da equação transcendental pelo procedimento do ponto médio. O valor que é obtido é θm=56.46º


Código: Seleccionar Todos
public class Ecuacion {
   static final double CERO=1e-10;
   static final double ERROR=0.001;
   static final int MAXITER=200;

public static void main(String[] args) {
   double aIni=50*Math.PI/180;
   double aFin=60*Math.PI/180;
   double raiz=puntoMedio(aIni, aFin);
   System.out.println(raiz*180/Math.PI);
}

static double puntoMedio(double a, double b) {
   double m, ym;
   int iter=0;
   do{
      m=(a+b)/2;
      ym=f(m);
      if(Math.abs(ym)<CERO) break;
      if(Math.abs((a-b)/m)<ERROR) break;

      if((f(a)*ym)<0) b=m;
      else a=m;
      iter++;
   }while(iter<MAXITER);
   if(iter==MAXITER){
      System.out.println("No se ha encontrado la raíz");
   }
   return m;
}

static double f(double x){
   double y=1.0-Math.sin(x)*Math.log((1+Math.sin(x))/Math.cos(x));
   return y;
}
}


Quando entrar na universidade aprenderei programar, infelizmente ainda falta 1 ano de colegio. /o/

MensagemEnviado: Sexta Jan 04, 2008 6:02 pm
por hexphreak
Ouch, até dá dó ver código Java :P

Zé, então usaste o mesmo método que eu :) Vou ver agora se completo a resolução analítica, tinha-me esquecido de um quadrado e por isso é que deu tudo mal :oops:

MensagemEnviado: Sexta Jan 04, 2008 7:15 pm
por jap
hexphreak Escreveu:Ouch, até dá dó ver código Java :P
(...)


Java já é mau, mas em espanhol (equacion, puntomedio??) é muito pior! :lol: :lol:

MensagemEnviado: Sexta Jan 04, 2008 7:27 pm
por jap
Aqui vai o equivalente código em Python (o código não é lá muito pitónico, para ser mais legível...)

Código: Seleccionar Todos
"""
The longest path

jap@2007 having fun for Quark!
"""

from math import *

def f(x):
    return sin(x)*log((1+sin(x))/cos(x))-1

def rad(x):
    return x*pi/180.

def deg(x):
    return x*180/pi
   
def bissect(f,a,b,eps=1E-12,maxit=500,show=0):
    """
    Finds a root of f(x) in the interval (a,b) to precision eps
    by sucessive bissection. Requires f(a)*f(b)<=0.
    """
    if f(a)*f(b) > 0:
        print "Error: no assurance of root in interval."
        return

    niter = 0
    while niter < maxit:
        if show: print niter,a,b
        mid = (a+b)/2.
        if f(mid)*f(a) > 0:
            a = mid
        else:
            b = mid
            if b-a < eps : return mid
        niter += 1
    else:
        print "Error: maximum nb. of iterations exceeded."


print deg(bissect(f,rad(10),rad(60),show=1))


Será dado um "aplauso" a quem conseguir resolver a equação pelo método IPF, em qualquer linguagem de programação :D

MensagemEnviado: Sexta Jan 04, 2008 8:20 pm
por hexphreak
jap Escreveu:Será dado um "aplauso" a quem conseguir resolver a equação pelo método IPF, em qualquer linguagem de programação :D

Percebe-se porquê :x O diabo dos negativos estraga tudo, mas eu chego lá...

MensagemEnviado: Sexta Jan 04, 2008 9:32 pm
por hexphreak
Viva a exponencial! :D

Código: Seleccionar Todos
#include <stdio>
#include <math>

int main(void)
{
    double sin_t = 0.8;
    int i;

    for(i = 0; i < 40; i++)
        sin_t = exp(2 / sin_t) * (1 - sin_t) - 1;

    printf("theta = %lf\n", asin(sin_t));
    return 0;
}

Código: Seleccionar Todos
$ gcc -lm -olong long.c
$ ./long
theta = -0.985515
$

Curiosamente, obtemos exactamente o simétrico do valor que se esperaria. Suponho que tenha a ver com a convergência da função que utilizei :)

Aqui está o raciocínio matemático: \sin \theta \ln{\frac{1 + \sin \theta}{\cos \theta}} = 1\Leftrightarrow \sin \theta \ln{\frac{1 + \sin \theta}{1 - \sin \theta}} = 2 \Leftrightarrow \sin \theta = e^{2 / \sin \theta} (1 - \sin \theta) - 1. Ou seja, manipulação trigonométrica, aplicar a exponencial e resolver em ordem ao seno :wink:

MensagemEnviado: Sexta Jan 04, 2008 10:57 pm
por jap
hexphreak Escreveu:Viva a exponencial! :D

Curiosamente, obtemos exactamente o simétrico do valor que se esperaria. Suponho que tenha a ver com a convergência da função que utilizei :)
(...)


:lol: Pois, mas nós não queremos o ângulo negativo, pois não? É solução matemática da equação, mas não física - o que é lançar um projéctil com ângulo negativo?

Consegues arranjar um IPF que convirja para o ângulo positivo? Ora tenta lá... :wink:

MensagemEnviado: Sexta Jan 04, 2008 11:19 pm
por hexphreak
\sin \theta = 1 - \frac{\sin \theta + 1}{e^{2/\sin \theta}}

Muito curioso :D

MensagemEnviado: Sexta Jan 04, 2008 11:32 pm
por jap
hexphreak Escreveu:
(...)
Aqui está o raciocínio matemático: \sin \theta\,\frac{1 + \sin \theta}{\cos \theta} = 1 \Leftrightarrow \sin \theta\,\ln{\frac{1 + \sin \theta}{1 - \sin \theta}} = 2 \Leftrightarrow \sin \theta = e^{2/\sin \theta} (1 - \sin \theta) - 1. Ou seja, manipulação trigonométrica, aplicar a exponencial e resolver em ordem ao seno :wink:

Desculpa, mas não percebi bem o teu raciocínio matemático. Verifica que não há um engano, algures... :?

MensagemEnviado: Sexta Jan 04, 2008 11:34 pm
por hexphreak
O que eu fiz foi resolver a equação que o Prof. postou em ordem ao seno, substituindo o co-seno e simplificando. Como a função não convergia no plano real (por causa do sinal do argumento para o logaritmo), apliquei a exponencial a ambos os lados da equação, e resolvi outra vez em ordem ao seno, mas parece que resolvi para o seno errado! :lol: Por isso experimentei com "outro seno" e funcionou, mas falta-me a base teórica.

Tinha-me esquecido de um logaritmo no post original :? Peço desculpa pela confusão.

MensagemEnviado: Sexta Jan 04, 2008 11:37 pm
por jap
Mas


\cos \theta = \sqrt{1- \sin^2 \theta}

Como é que a partir daí efectuaste simplificações para chegar à tua expressão final? :roll:

MensagemEnviado: Sexta Jan 04, 2008 11:48 pm
por hexphreak
\sin \theta = 1/\ln{\frac{1+\sin \theta}{\sqrt{1-\sin^2 \theta}}}

Ora, 1 - \sin^2 \theta = (1 - \sin \theta)(1 + \sin \theta), logo:

\sqrt{\frac{(1+\sin \theta)^2}{(1-\sin \theta)(1+\sin \theta)}} = \sqrt{\frac{1+\sin \theta}{1-\sin \theta}}

E podemos agora substituir na expressão de cima:

\sin \theta = 1/\ln{\sqrt{\frac{1+\sin \theta}{1-\sin \theta}}} = 2/\ln{\frac{1+\sin \theta}{1-\sin \theta}}

Logo

(\frac{1+\sin \theta}{1-\sin \theta})^{\sin \theta} = e^2 \Leftrightarrow \sin \theta = 1 - \frac{1+\sin \theta}{e^{2/\sin \theta}}

MensagemEnviado: Sexta Jan 04, 2008 11:55 pm
por jap
OK, desculpa, fiquei inicialmente confundido porque faltava um logaritmo algures na tua expressão e eu não reparei na gralha :oops:
Está tudo claro, e as substituições que fizeste para dar a volta à expressão são fixes. Obrigado! Este novo IPF funciona, ou seja converge para a raiz correcta? Ora mostra aí à malta, começando de vários pontos de partida (razoáveis!) :D

MensagemEnviado: Sábado Jan 05, 2008 12:14 am
por hexphreak
Funciona pois :D Modifiquei ligeiramente o programa para mostrar as iterações e para correr vários valores iniciais para o seno, podem encontrar os resultados aqui (em radianos) :)

MensagemEnviado: Sábado Jan 05, 2008 12:31 am
por jap
ou aqui :lol:

Código: Seleccionar Todos
"""
The longest path
IPF method

jap@2007 having fun for Quark!
"""

from math import *

def g(x):
    return 1-(x+1)/exp(2/x)

def rad(x):
    return x*pi/180.

def deg(x):
    return x*180/pi
   
angles = [10*n for n in range(1,10)]
print "Trial angles=",angles

for angle in angles:
    x = sin(rad(angle))
    for n in range(50):
        x = g(x)
    print "Initial angle",angle,"converges to",deg(asin(x))

print "Big hand of applause for hexphreak that figured out a very stable IPF! :)"


Output:

Código: Seleccionar Todos
>>>
>>>
Trial angles= [10, 20, 30, 40, 50, 60, 70, 80, 90]
Initial angle 10 converges to 56.4658351274
Initial angle 20 converges to 56.4658351274
Initial angle 30 converges to 56.4658351274
Initial angle 40 converges to 56.4658351274
Initial angle 50 converges to 56.4658351274
Initial angle 60 converges to 56.4658351275
Initial angle 70 converges to 56.4658351275
Initial angle 80 converges to 56.4658351275
Initial angle 90 converges to 56.4658351275
Big hand of applause for hexphreak that figured out a very stable IPF! :)
>>> >>>


:hands: