Quiet@Remove["Global`*"];
$HistoryLength=2;
SetDirectory[NotebookDirectory[]];

?_

NestList[f,x0,5]

(*the simpliest iteration*)
Quiet@Remove["Global`*"];
(*we solve eqn Cos[x]==x*)
xNext[xPrevious_]:=Cos[xPrevious];
xStart=0.123;
sol=NestList[xNext,xStart,20]
ListPlot[sol]

(((-3))^3)

((-27)^(1/3))//N

Quiet@Remove["Global`*"];
(*we solve eqn 5(x^3)+2(x^2)+x+3==0*)
xNext[xPrevious_]:=-(3.-2(xPrevious^2)-5(xPrevious^3));
xStart=0;
sol=NestList[xNext,xStart,10]

Solve[5(x^3)+2(x^2)+x+3==0]//N

FindRoot[5(x^3)+2(x^2)+x+3==0,{x,0.2+0.7I}]

Quiet@Remove["Global`*"];
(*we solve eqn 5(x^3)+2(x^2)+x+3==0*)
xNext[xPrevious_]:=(((1./5)*(-2(xPrevious^2)-xPrevious-3))^(1/3));
xStart=-5;
sol=NestList[xNext,xStart,10]

Quiet@Remove["Global`*"];
(*Newton Method*)
(*we solve eqn f[x]==0*)
f[x_]:=0.5x+0.1Sin[x];
xNext[xPrevious_]:=xPrevious-((f[xPrevious])/(f'[xPrevious]));
xStart=-5;
sol=NestList[xNext,xStart,5]

?/;

f[-6]*f[6]

Quiet@Remove["Global`*"];
(*bisection*)
f[x_]:=0.5x+0.1Sin[x]+1;
nextInterval[{a_,b_}]:={a,((a+b)/2)}/;f[a]*f[((a+b)/2)]<0;
nextInterval[{a_,b_}]:={((a+b)/2),b}/;f[b]*f[((a+b)/2)]<0;
nextInterval[{a_,b_}]:={a,a}/;f[a]==0;
nextInterval[{a_,b_}]:={b,b}/;f[b]==0;
startingInterval={-6.,6.};
NestList[nextInterval,startingInterval,10]

(*optimization*)
ConstrainedMin[2 x-3 y,{x+y<10,x-y>2,x>1},{x,y}]

(*for linear problems use LinearPrograming*)
LinearProgramming[{2,-3},(({{-1, -1}, {1, -1}, {1, 0}})),{-10,2,1}]

(x^2)+(y^3)/.x->Pi

?/.

(*patterns and functions*)
Remove["Global`*"]
f[x:_]:=(x^2);
f[5]

Remove["Global`*"]
f[x_]:=(x^2);
f[5]

(*for other problems we use NMinimize , NMaximize, FindMinimum, FindMaximum*)
(*example of constrained minimum*)
(*clasic solution, we nearly newer use it*)
expr=x+(y^2)+Sin[x+y];
constr=(x^2)+(y^2)==1;
g=constr/.lhs_==rhs_:>lhs-rhs;
lagrange=expr+λ*g;
eqns={D[lagrange,x]==0,D[lagrange,y]==0,D[lagrange,λ]==0}
FindRoot[eqns,{x,0.5},{y,0.5},{λ,0.5}]
FindRoot[eqns,{x,-0.5},{y,-0.5},{λ,1}]

expr=x+(y^2)+Sin[x+y];
constr=(x^2)+(y^2)==1;
NMinimize[{expr,constr},{x,y}]
NMaximize[{expr,constr},{x,y}]

Plot[Tanh[20(x-2)]-Tanh[20(x-4)],{x,-3,6}]

Remove["Global`*"]
expr=x+(y^2)+Sin[x+y];
constr1=2<x<4;
constr2=5<y<6;
NMinimize[{expr,constr1,constr2},{x,y}]

penaltyY=(((y-5.5))^2)*(2-(Tanh[200(y-5)]-Tanh[200(y-6)]))
Plot[penaltyY,{y,4,7}]

penaltyX=(((x-3))^2)*(2-(Tanh[200(x-2)]-Tanh[200(x-4)]))
Plot[penaltyX,{x,-3,9}];

newExpr=expr+penaltyX+penaltyY;
FindMinimum[newExpr,{{x,0},{y,0}}]

Head[5]

Head[5.]

(7/13)-(11./19)

Quiet@Remove["Global`*"];
(*optimization for solution of equations*)
eqns={x+y==1,x-2y==5.};
Solve[eqns]

forMin=Total[eqns/.lhs_==rhs_:>(((lhs-rhs))^2)]
NMinimize[forMin,{x,y}]

{a,b,c}=={d,e,f}
Thread[{a,b,c}=={d,e,f}]

matrix.{x,y,z}==vector

Quiet@Remove["Global`*"];
(*overdetermined systems*)
matrix=Table[Random[],{5},{3}];
vector=Table[Random[],{5}];
eqns=Thread[matrix.{x,y,z}==vector];
forMin=Total[eqns/.lhs_==rhs_:>(((lhs-rhs))^2)];
NMinimize[forMin,{x,y,z}]
b=Transpose[matrix].vector;
A=Transpose[matrix].matrix;
LinearSolve[A,b]//N

{MatrixForm[matrix],MatrixForm[vector]}

matrix.{x,y,z}==vector
Thread[matrix.{x,y,z}==vector]

LinearSolve

MatrixForm[Transpose[matrix].matrix]
Det[Transpose[matrix].matrix]//N

Quiet@Remove["Global`*"];
(*overdetermined systems*)
matrix=Table[i+j+(Sqrt[i*j]),{i,5},{j,3}];
vector=Table[Prime[i],{i,5}];
{MatrixForm[matrix],MatrixForm[vector]};
eqns=Thread[matrix.{x,y,z}==vector];
forMin=Total[eqns/.lhs_==rhs_:>(((lhs-rhs))^2)];
NMinimize[forMin,{x,y,z}]
b=Transpose[matrix].vector;
A=Transpose[matrix].matrix;
LinearSolve[A,b]//N

Quiet@Remove["Global`*"];
points={a,b,c,d}={{-1,-1.5},{3,1.},{1.,-2},{-1.5,-1.5}};
myPoint={myX,myY}={5,5};
distance[point:{x_,y_}]:=(Sqrt[(((myX-x))^2)+(((myY-y))^2)]);
distances=distance/@points

distance2[point:{x_,y_}]:=(Sqrt[(((myx-x))^2)+(((myy-y))^2)]);
error[point_,distance_]:=(((distance2[point]-distance))^2);
totalerror=error[a,8.9]+error[b,4.5]+error[c,8.1]+error[d,9.2];
NMinimize[totalerror,{myx,myy}]

Quiet@Remove["Global`*"];
(*a quadratic fit*)
originalFunction[x_]:=(x^2)-2x+1;
SeedRandom[1];
data=Table[{x,originalFunction[x]+RandomReal[{-0.5,0.5}]},{x,-2,3,0.2}];
pl1=ListPlot[data,PlotStyle->{Red,PointSize[0.015]}]

(*least squares*)

pl2=Plot[originalFunction[x],{x,-2,3}];
Show[pl1,pl2];
error[{x_,y_}]:=(((y-(a*(x^2)+b*x+c)))^2);
totalError=(∑_(i=1)^(Length[data]))error[data[[i]]];
NMinimize[totalError,{a,b,c}]

totalError//Simplify

eqns={D[totalError,a]==0,D[totalError,b]==0,D[totalError,c]==0}//Simplify
Solve[eqns]

pl1=ListPlot[data,PlotStyle->{Red,PointSize[0.015]}];
fit=Fit[data,{1,x,(x^2),(x^3),(x^4),(x^5),(x^6),(x^7)},x]
pl3=Plot[fit,{x,-2,3}];
Show[pl1,pl3]

Quiet@Remove["Global`*"];
(*example Fourier Series*)
f=50;
T=(1/f);
u[t_]:=Tanh[2*Sin[2*Pi*f*(t-0.123T)]];
pl1=Plot[u[t],{t,0,2T}];
data=Table[{t,u[t]},{t,0,T,(T/300.)}];
n=7;
fces=Table[{Sin[i*2*Pi*f*t],Cos[i*2*Pi*f*t]},{i,1,n,1}]//Flatten;
fit=Fit[data,fces,t];
pl2=Plot[fit,{t,0,2T},PlotStyle->{Red,Dashed}];
Show[pl1,pl2];

Quiet@Remove["Global`*"];
(*nonlinear problems*)
f[t_]:=3-2*(E^((-t)/3));
pl1=Plot[f[t],{t,0,10},AxesOrigin->{0,0}];
data=Table[{t,f[t]},{t,0,10,0.2}];
err[{t_,u_}]:=(((u-(a+b*(E^((-t)/τ)))))^2);
totalError=Total[err/@data];
NMinimize[totalError,{a,b,τ}];
params=FindFit[data,a+b*(E^((-t)/τ)),{a,b,τ},t];
fit=a+b*(E^((-t)/τ))/.params;
pl2=ListPlot[data,PlotStyle->{Red,PointSize[0.015]}];
pl3=Plot[fit,{t,0,10}];
Show[pl2,pl3];

Quiet@Remove["Global`*"];
(*numerical integration*)
(*Simpson's rule, derivation*)
expr[t_]:=a*(t^2)+b*t+c;
eqns={expr[t0-Δt]==fMinus,expr[t0+Δt]==fPlus,expr[t0]==f0};
fit=expr[t]/.Solve[eqns,{a,b,c}][[1]];
integral=Integrate[fit,{t,t0-Δt,t0+Δt}]//Simplify

{OddQ[5],OddQ[4]}
{EvenQ[5],EvenQ[4]}

Integrate[Sin[t],{t,0,Pi}]

Quiet@Remove["Global`*"];
(*numerical integration*)
(*Simpson's rule, example*)
(*nPoints must be odd integer*)
simpson[f_,tmin_,tmax_,nPoints_?OddQ]:=Module[{tValues,Δt,functionValues,triples},
Δt=((tmax-tmin)/(nPoints-1));
tValues=Range[tmin,tmax,Δt];
functionValues=f/@tValues;
triples=Partition[functionValues,3,2];
((Δt)/3)Total[triples/.{a_?NumericQ,b_,c_}:>a+4b+c]
];

((simpson[Sin,0.,Pi,5]-2)/2)*100

list=f[#]&/@Range[11]
triples=Partition[list,3,2]
Total[triples/.{a_,b_,c_}:>a+4b+c]

f/@{1,2,3,4}

Quiet@Remove["Global`*"];
(*numerical integration*)
(*Simpson's rule, example*)
(*nPoints must be odd integer*)
simpson[f_,tmin_,tmax_,nIntervals_?EvenQ]:=Module[{tValues,Δt,functionValues,triples},
Δt=((tmax-tmin)/nIntervals);
tValues=Range[tmin,tmax,Δt];
functionValues=f/@tValues;
triples=Partition[functionValues,3,2];
((Δt)/3)Total[triples/.{a_?NumericQ,b_,c_}:>a+4b+c]
];

((simpson[Sin,0.,Pi,4]-2)/2)*100

Quiet@Remove["Global`*"];
(*numerical integration*)
(*Simpson's rule, example*)
(*nPoints must be odd integer*)
simpson[f_,tmin_,tmax_,nPoints_?OddQ]:=Module[{tValues,Δt,functionValues,triples},
Δt=((tmax-tmin)/(nPoints-1));
tValues=Range[tmin,tmax,Δt];
functionValues=f/@tValues;
triples=Partition[functionValues,3,2];
((Δt)/3)Total[triples/.{a_?NumericQ,b_,c_}:>a+4b+c]
];

error[nPoints_?OddQ]:=((simpson[Sin,0.,Pi,nPoints]-2)/2)*100;
nMin=7;
nMax=51;
data=Table[{n,error[n]},{n,nMin,nMax,2}];
pl1=ListPlot[data,PlotRange->All];
fit=Fit[data,{(1/(n^4))},n];
pl2=Plot[fit,{n,nMin,nMax},PlotRange->All];
Show[pl1,pl2,PlotLabel->fit]

Quiet@Remove["Global`*"];
(*ODE, LaplaceTransform*)
LaplaceTransform[y'[t],t,s]/.{LaplaceTransform[y[t],t,s]->Y[s]}

LaplaceTransform[y''[t],t,s]/.{LaplaceTransform[y[t],t,s]->Y[s]}

eqns={s*x1==-2x1+x2+u,s*x2==x1-3x2};
Solve[eqns,{x1,x2}]

Quiet@Remove["Global`*"];
noise:=0.001RandomReal[{-1,1}];
data=Table[{x,noise+Sin[x]},{x,0,4Pi,0.02}];
signal=Interpolation[data];
{Plot[signal[x],{x,0,4Pi}],Plot[signal'[x],{x,0,4Pi}],Plot[signal''[x],{x,0,4Pi}]}

NestList[f,x0,5]

Quiet@Remove["Global`*"];
(*Euler method*)
derivative[{t_,y_}]:=-y+Sin[t];
tStart=0.;
yStart=10;
tEnd=20;
nSteps=10;
Δt=((tEnd-tStart)/nSteps);
step[point:{t_,y_}]:={t+Δt,y+derivative[point]*Δt};
sol=NestList[step,{tStart,yStart},nSteps];
ListPlot[sol,PlotRange->All]

derivative[{t_,y_}]:=-y+Sin[t];
tStart=0.;
yStart=10;
tEnd=20;
nSteps=10;
Δt=((tEnd-tStart)/nSteps);
RungeKuttaStep[{t_,y_}]:=Module[
{k1,k2,k3,k4},
k1=derivative[{t,y}];
k2=derivative[{t+((Δt)/2),y+((Δt)/2)*k1}];
k3=derivative[{t+((Δt)/2),y+((Δt)/2)*k2}];
k4=derivative[{t+Δt,y+Δt*k3}];
{t+Δt,y+((Δt)/6)*(k1+2k2+2k3+k4)}
];
sol=NestList[RungeKuttaStep,{tStart,yStart},nSteps];
ListPlot[sol,PlotRange->All]

RungeKuttaStep[{t_,y_}]:=Module[
{k1,k2,k3,k4},
k1=derivative[{t,y}];
k2=derivative[{t+((Δt)/2),y+((Δt)/2)*k1}];
k3=derivative[{t+((Δt)/2),y+((Δt)/2)*k2}];
k4=derivative[{t+Δt,y+Δt*k3}];
{t+Δt,y+((Δt)/6)*(k1+2k2+2k3+k4)}
];

Solve[x==1]

Solve[(x^2)==1]

DSolve[((y'[t])^2)==-y[t],y[t],t]

Quiet@Remove["Global`*"];
(*original eqn*)
eqn={y''[t]+2y'[t]+y[t]==10*Sin[10t],y[0]==0,y'[0]==1};
analSol=DSolve[eqn,y[t],t][[1]];
tmax=10;
pl1=Plot[y[t]/.analSol,{t,0,tmax}];
numSol=NDSolve[eqn,y[t],{t,0,tmax}][[1]];
Plot[y[t]/.numSol,{t,0,tmax}];
derivative[{t_,{x1_,x2_}}]:={x2,10*Sin[10t]-x1-2x2};

nSteps=50;
Δt=((tmax-0)/nSteps);
start={0.,{0.,1}};
RungeKuttaStep[{t_,y_}]:=Module[
{k1,k2,k3,k4},
k1=derivative[{t,y}];
k2=derivative[{t+((Δt)/2),y+((Δt)/2)*k1}];
k3=derivative[{t+((Δt)/2),y+((Δt)/2)*k2}];
k4=derivative[{t+Δt,y+Δt*k3}];
{t+Δt,y+((Δt)/6)*(k1+2k2+2k3+k4)}
];
sol=NestList[RungeKuttaStep,start,nSteps];
data=sol/.{t_,{x1_,x2_}}:>{t,x1};
Show[pl1,ListPlot[data,PlotStyle->Red]]

(x^2)/.{x->2}

Plot[Sin[x],{x,0,6}]

Plot[{Sin[x]},{x,0,6}]

(x^2)/.{{x->2}}

Quiet@Remove["Global`*"];
(*original eqn*)
eqn={y''[t]+2y'[t]+y[t]==10*Sin[10t],y[0]==0,y'[0]==1};
analSol=DSolve[eqn,y[t],t][[1]];
tmax=10;
pl1=Plot[y[t]/.analSol,{t,0,tmax}];
numSol=NDSolve[eqn,y[t],{t,0,tmax}][[1]]
solution[t_]=y[t]/.numSol;
solution[0.5tmax]

Quiet@Remove["Global`*"];
$HistoryLength=2;
SetDirectory[NotebookDirectory[]];

(*shooting method*)
(*example, heat transfer in a wall with nonconstant thermal conductivity*)
λ[T_]:=0.4+0.01*T;
eqn=D[λ[T[x]]*T'[x],x]==0;
T0=-5;
d=0.45;
Td=20;
sol[derivative_]:=NDSolve[{eqn,T[0]==T0,T'[0]==derivative},T,{x,0,d}][[1]];
TxEquald[derivative_?NumericQ]:=T[d]/.sol[derivative];
solder=FindRoot[TxEquald[der]==Td,{der,10}]
solution=sol[der/.solder];
Plot[T[x]/.solution,{x,0,d}]

Quiet@Remove["Global`*"];
(*numerical derivatives*)
f[x_]:=x*Sin[x+1];
numDer[x_,Δx_]:=((f[x+Δx]-f[x])/(Δx));
numDer2[x_,Δx_]:=((f[x+Δx]-f[x-Δx])/(2Δx));
Plot[{f'[x]-numDer[x,0.1],f'[x]-numDer2[x,0.1]},{x,0,2Pi},PlotStyle->{Red,Blue}]

Quiet@Remove["Global`*"];
(*numerical derivatives*)
f[x_]:=x*Sin[x+1];
Δx=0.1;
numDer[x_]:=((f[x+Δx]-f[x])/(Δx));
num2ndDer[x_]:=((numDer[x+Δx]-numDer[x])/(Δx));
numDer2[x_]:=((f[x+Δx]-f[x-Δx])/(2Δx));
num2ndDer2[x_]:=((numDer2[x+Δx]-numDer2[x-Δx])/(2Δx));
num2Der[x_]:=((f[x+Δx]+f[x-Δx]-2*f[x])/(((Δx))^2));

Plot[{f''[x]-num2Der[x],f''[x]-num2ndDer2[x]},{x,0,2Pi},PlotStyle->{Red,Blue}]

(*dangers in ODEs*)
Quiet@Remove["Global`*"];
lambda=Module[{eqns},
eqns={a-b==0.04,a+b==0.8};
a+b*Tanh[((x-0.1)/0.00001)]/.Solve[eqns][[1]]
];
Plot[lambda,{x,0,0.55}]

lambda

λ=lambda;
eqn=eqn=D[λ*T'[x],x]==0;
T0=-5;
d=0.55;
Td=20;

λ=lambda;
eqn=eqn=D[λ*T'[x],x]==0;
T0=-5;
d=0.55;
Td=20;
sol[derivative_]:=NDSolve[{eqn,T[0]==T0,T'[0]==derivative},T,{x,0,d},MaxStepSize->0.0001][[1]];
TxEquald[derivative_?NumericQ]:=T[d]/.sol[derivative];
solder=Quiet@FindRoot[TxEquald[der]==Td,{der,100}];
solution=sol[der/.solder];
Plot[T[x]/.solution,{x,0,d},PlotRange->All]

derλ=D[λ,x];
Plot[derλ,{x,0,0.55},PlotRange->All,PlotPoints->300]

Quiet@Remove["Global`*"];
SetDirectory[NotebookDirectory[]];
(*frequency estimation*)
f=49.2;
T=(1/f);
ω=2*Pi*f;
u[t_]:=230*(Sqrt[2])*Sin[ω*t+0.123]+10*Sin[3*ω*t+0.523];
data=Table[{t,u[t]},{t,0.,5T,0.00713T}];
ListPlot[data]
Export["dataFrequency.csv",data]

list=(Subscript[a, #])&/@Range[10]
Partition[list,3,1]

Partition[Take[data,5],2,1]

data=Import["dataFrequency.csv"];
plData=ListPlot[data];
pairsOfPoints=Partition[data,2,1];
test[{{t1_,u1_},{t2_,u2_}}]:=u1*u2<0;
zeroCrossingPairs=Select[pairsOfPoints,test];
plSelected=ListPlot[zeroCrossingPairs,PlotStyle->Red];
Show[plData,plSelected];
timeOfCrossing[{{t1_,u1_},{t2_,u2_}}]:=Module[{fit},
fit=Fit[{{t1,u1},{t2,u2}},{1,t},t];
t/.Solve[fit==0][[1]]
];
timeOfCrossing[{{1,2},{2,-0.5}}]

g/@{a,b,c}

list=(Subscript[t, #])&/@Range[10]
Partition[list,3,1]
%/.{t1_,t2_,t3_}:>t3-t1

timeOfCrossing[{{t1_,u1_},{t2_,u2_}}]:=Module[{fit},
fit=Fit[{{t1,u1},{t2,u2}},{1,t},t];
Solve[fit==0][[1,1,2]]
];
times=timeOfCrossing/@zeroCrossingPairs;
periods=Partition[times,3,1]/.{t1_,t2_,t3_}:>t3-t1;
period=Mean[periods];
Print["frequency= ",(1/period)," Hz"]