Got Questions? Get Answers.
Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
optimization of fiber bragg grating using nelder-mead simplex algoritm

Subject: optimization of fiber bragg grating using nelder-mead simplex algoritm

From: kuffal al

Date: 4 Apr, 2010 14:54:05

Message: 1 of 5

can someone help me to optimize the fiber bragg grating(fbg) using nelder mead simplex algoritm. i have both codes (both codes working well) the fbg n nelder mead but i fail to integrate them. the parameter that i want to optimize are the length (L), period (A) and the induced index change (dn). so here are the codes

(1)the fbg code
%==========================================================================
% matlab script for solving the coupled-mode equation using
% the transfer matrix method
%==========================================================================

clear all
close all
clc

%==========================================================================
% fibre simulation parameters
%==========================================================================

L = 0.0100; % Grating Length

dn = 1.0e-4; % Induced index change

design_lambda = 1550e-9; % design wavelength

%==========================================================================
% finding neff(effective reflective index) for the designed fibre
%==========================================================================

format long
n=40;
a=0.5;
increment=2;
tolerance=1e-4;

nco=1.447; % nco-core index
ncl=1.444; % ncl-cladding index
xstart=((ncl-0.0112)+nco)/2;
x=xstart;
dx=increment;

%==========================================================================
% to find the root of neff(effective reflective index)
%==========================================================================


for m=1:n
    neff= sqrt((ncl^2)+(a*x)*(nco^2-ncl^2));
    while dx/x>tolerance
        if neff~=sqrt((ncl^2)+(a*(x+dx))*(nco^2-ncl^2));
            dx=dx/2;
       else
            x=x+dx;
        end
    end
    route(m)=x;
    dx=increment;
    x=1.0001*x;
end

neff=route(40:end);

V =(2*pi/design_lambda)*a*sqrt(nco^2-ncl^2); % normalize frequency

A=design_lambda/(2*neff); % Grating period 5.36e-7;

lambda=(1549.50e-9:0.005e-9:1550.50e-9); % Spectra range

 
%==========================================================================
% Coupled-mode Theory
%==========================================================================

v=1; % visibility, v assumed to be 1

k=(pi./lambda)*v*dn; % AC coupling coefficient

% DC(period averaged) coupling coefficient

% for stronger grating, sigma_dc assumed to be 0

sigma_dc =0; %(2*pi./lambda).*dn;

 % general DC coupling coefficient

sigma = 2*pi*neff.*(1./lambda - 1/(2*neff*A)) + sigma_dc;

 
%==========================================================================
% Power Reflection coefficients(R)
%==========================================================================

K2 = k.*k;

sigma2 = sigma.*sigma;

difference = K2-sigma2;

up = sinh(L*sqrt(difference)).^2;

down = cosh(L*sqrt(difference)).^2 - sigma2./K2;

R = up./down;

Rt = transpose(R);

%==========================================================================
% saving the points
%==========================================================================

save 'fileToRead1.txt' Rt -ascii

%==========================================================================
% plots for reflection spectra for Bragg reflector
%==========================================================================

figure(1)

plot(lambda,R,'b');

axis( [1549.50e-9, 1550.50e-9, 0, 1]);

title('Reflection Spectra of Bragg grating');

xlabel('Wavelength(m)');

ylabel('Reflectivity(p.u)');

figure(2)
stairs(lambda,R,'r');figure(gcf)
axis( [1549.50e-9, 1550.50e-9, 0, 1]);

--------------------------------------------------------------------------------------------------------------------

(2) the nelder-mead code

%*******************************************
%* Nelder and Mead Technique *
%* Siamak Faridani *
%* University of Oklahoma *
%* 16 Oct 2005 *
%* www.PitchUp.com *
%*******************************************


clc;
clear all;
close all;
N=2;
alpha=2;
beta=.25;
gamma=2.5;
Precision=.00001

%Insert your starting point and the function here
x01=1.5505*1e-6;
x02=1.550e-6;
func=@Modified_fbg


%do not modify the following lines
%==============================

x=[[x01,x02]];



delta1=(((N+1)^.5+N-1)/(N*sqrt(2)))*alpha;

delta2=(((N+1)^.5-1)/(N*sqrt(2)))*alpha
x=[x;[x(1,1)+delta1,x(1,2)+delta2];[x(1,1)+delta2,x(1,2)+delta1]];
x_archive=[x [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]];
l=0
x
disp('start================');
Precision_inverse=Precision^-1;
while ~(isequal(round(x(1,:)*Precision_inverse),round(x(2,:)*Precision_inverse),round(x(3,:)*Precision_inverse)))
    l=l+1;
    values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
   
    sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
    xc=(1/N)*(x(sorted(1,2),:)+x(sorted(1,3),:));
    disp('trying Normal Reflection');
    xnew=x(sorted(1,1),:)+(1+alpha)*(xc-x(sorted(1,1),:))
    isitok=0;
    func(xnew(1,1),xnew(1,2))
    while isitok==0
        
    if ((func(xnew(1,1),xnew(1,2))>values(sorted(1,3),:)) && (func(xnew(1,1),xnew(1,2))<values(sorted(1,2),:)))
        disp('f(l)<f(x_new)<f(g)')
        isitok=1;
    elseif (func(xnew(1,1),xnew(1,2))<values(sorted(1,3),:))
        disp('lets Expand...');
        xnew1=x(sorted(1,1),:)+(1+gamma)*(xc-x(sorted(1,1),:))
        func(xnew1(1,1),xnew1(1,2))
        if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
            xnew=xnew1
        end
                isitok=1;
    elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) && (func(xnew(1,1),xnew(1,2))>values(sorted(1,1),:))
        disp('lets Contract with teta=-beta...');
        xnew1=x(sorted(1,1),:)+(1-beta)*(xc-x(sorted(1,1),:))
        func(xnew1(1,1),xnew1(1,2))
        if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
            xnew=xnew1
        end
                isitok=1;
    elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) && (func(xnew(1,1),xnew(1,2))<values(sorted(1,1),:))
        disp('lets Contract with teta=beta...');
        xnew1=x(sorted(1,1),:)+(1+beta)*(xc-x(sorted(1,1),:))
        func(xnew1(1,1),xnew1(1,2))
        if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
            xnew=xnew1
        end
        isitok=1;
    end
    end
    x(sorted(1,1),:)=xnew;
    x
    disp('Finished================');
    x_archive=[x_archive;[x [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]]];
end

subplot(1,2,1)
plot(x_archive(:,1),x_archive(:,2),'o-r')
subplot(1,2,2)
patch(x_archive(:,1),x_archive(:,2),x_archive(:,3),'r')
view(3)

values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
   
    sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
disp(sprintf('Number of iteration= %10d\n',l))
disp(sprintf('minimum of f(x) is %10.5f that occurs at x=[%10.5f,%10.5f]\n',min(values),x(find(values==min(values)),1),x(find(values==min(values)),2)))

Subject: optimization of fiber bragg grating using nelder-mead simplex

From: Alan Weiss

Date: 4 Apr, 2010 16:22:27

Message: 2 of 5

You might want to use some built-in MATLAB functions.
fminsearch performs the Nelder-Mead simplex algorithm.
It accepts a single vector of inputs, say x = [x(1),x(2),x(3)], and
works on your objective function that should take x and return f(x),
where f(x) is the function to minimize. For your problem, take x =
[L,A,dn], and write your objective function in terms of x.

It seems to me that your code is doing some very inefficient
root-finding as well. Take a look at the fzero function for a solid
root-finding implementation.

fminsearch reference:
http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fminsearch.html

fminsearch example:
http://www.mathworks.com/access/helpdesk/help/techdoc/math/bsgpq6p-10.html

fzero reference:
http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fzero.html

fzero example:
http://www.mathworks.com/access/helpdesk/help/techdoc/math/bsgpq6q-44.html

Alan Weiss
MATLAB mathematical toolbox documentation

kuffal al wrote:
> can someone help me to optimize the fiber bragg grating(fbg) using
> nelder mead simplex algoritm. i have both codes (both codes working
> well) the fbg n nelder mead but i fail to integrate them. the parameter
> that i want to optimize are the length (L), period (A) and the induced
> index change (dn). so here are the codes
>
> (1)the fbg code
> %==========================================================================
> % matlab script for solving the coupled-mode equation using % the
> transfer matrix method
> %==========================================================================
>
> clear all
> close all
> clc
>
> %==========================================================================
> % fibre simulation parameters
> %==========================================================================
>
> L = 0.0100; % Grating Length
>
> dn = 1.0e-4; % Induced index change
>
> design_lambda = 1550e-9; % design wavelength
>
> %==========================================================================
> % finding neff(effective reflective index) for the designed fibre
> %==========================================================================
>
> format long
> n=40;
> a=0.5;
> increment=2;
> tolerance=1e-4;
>
> nco=1.447; % nco-core index
> ncl=1.444; % ncl-cladding index
> xstart=((ncl-0.0112)+nco)/2;
> x=xstart;
> dx=increment;
>
> %==========================================================================
> % to find the root of neff(effective reflective index)
> %==========================================================================
>
>
> for m=1:n
> neff= sqrt((ncl^2)+(a*x)*(nco^2-ncl^2));
> while dx/x>tolerance
> if neff~=sqrt((ncl^2)+(a*(x+dx))*(nco^2-ncl^2));
> dx=dx/2;
> else
> x=x+dx;
> end
> end
> route(m)=x;
> dx=increment;
> x=1.0001*x;
> end
> neff=route(40:end);
>
> V =(2*pi/design_lambda)*a*sqrt(nco^2-ncl^2); % normalize frequency
>
> A=design_lambda/(2*neff); % Grating period 5.36e-7;
>
> lambda=(1549.50e-9:0.005e-9:1550.50e-9); % Spectra range
>
>
> %==========================================================================
> % Coupled-mode Theory
> %==========================================================================
>
> v=1; % visibility, v assumed to be 1
>
> k=(pi./lambda)*v*dn; % AC coupling coefficient
>
> % DC(period averaged) coupling coefficient
>
> % for stronger grating, sigma_dc assumed to be 0
>
> sigma_dc =0; %(2*pi./lambda).*dn;
>
> % general DC coupling coefficient
>
> sigma = 2*pi*neff.*(1./lambda - 1/(2*neff*A)) + sigma_dc;
>
>
> %==========================================================================
> % Power Reflection coefficients(R)
> %==========================================================================
>
> K2 = k.*k;
>
> sigma2 = sigma.*sigma;
>
> difference = K2-sigma2;
>
> up = sinh(L*sqrt(difference)).^2;
>
> down = cosh(L*sqrt(difference)).^2 - sigma2./K2;
>
> R = up./down;
>
> Rt = transpose(R);
>
> %==========================================================================
> % saving the points
> %==========================================================================
>
> save 'fileToRead1.txt' Rt -ascii
>
> %==========================================================================
> % plots for reflection spectra for Bragg reflector
> %==========================================================================
>
> figure(1)
>
> plot(lambda,R,'b');
>
> axis( [1549.50e-9, 1550.50e-9, 0, 1]);
>
> title('Reflection Spectra of Bragg grating');
>
> xlabel('Wavelength(m)');
>
> ylabel('Reflectivity(p.u)');
>
> figure(2)
> stairs(lambda,R,'r');figure(gcf)
> axis( [1549.50e-9, 1550.50e-9, 0, 1]);
>
> --------------------------------------------------------------------------------------------------------------------
>
>
> (2) the nelder-mead code
>
> %*******************************************
> %* Nelder and Mead Technique *
> %* Siamak Faridani *
> %* University of Oklahoma *
> %* 16 Oct 2005 *
> %* www.PitchUp.com *
> %*******************************************
>
>
> clc;
> clear all; close all;
> N=2;
> alpha=2;
> beta=.25;
> gamma=2.5;
> Precision=.00001
>
> %Insert your starting point and the function here
> x01=1.5505*1e-6;
> x02=1.550e-6;
> func=@Modified_fbg
>
>
> %do not modify the following lines
> %==============================
>
> x=[[x01,x02]];
>
>
>
> delta1=(((N+1)^.5+N-1)/(N*sqrt(2)))*alpha;
>
> delta2=(((N+1)^.5-1)/(N*sqrt(2)))*alpha
> x=[x;[x(1,1)+delta1,x(1,2)+delta2];[x(1,1)+delta2,x(1,2)+delta1]];
> x_archive=[x
> [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]];
> l=0
> x
> disp('start================');
> Precision_inverse=Precision^-1;
> while
> ~(isequal(round(x(1,:)*Precision_inverse),round(x(2,:)*Precision_inverse),round(x(3,:)*Precision_inverse)))
>
> l=l+1;
> values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
>
> sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
>
> xc=(1/N)*(x(sorted(1,2),:)+x(sorted(1,3),:));
> disp('trying Normal Reflection');
> xnew=x(sorted(1,1),:)+(1+alpha)*(xc-x(sorted(1,1),:))
> isitok=0;
> func(xnew(1,1),xnew(1,2))
> while isitok==0
> if ((func(xnew(1,1),xnew(1,2))>values(sorted(1,3),:)) &&
> (func(xnew(1,1),xnew(1,2))<values(sorted(1,2),:)))
> disp('f(l)<f(x_new)<f(g)')
> isitok=1;
> elseif (func(xnew(1,1),xnew(1,2))<values(sorted(1,3),:))
> disp('lets Expand...');
> xnew1=x(sorted(1,1),:)+(1+gamma)*(xc-x(sorted(1,1),:))
> func(xnew1(1,1),xnew1(1,2))
> if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
> xnew=xnew1
> end
> isitok=1;
> elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) &&
> (func(xnew(1,1),xnew(1,2))>values(sorted(1,1),:))
> disp('lets Contract with teta=-beta...');
> xnew1=x(sorted(1,1),:)+(1-beta)*(xc-x(sorted(1,1),:))
> func(xnew1(1,1),xnew1(1,2))
> if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
> xnew=xnew1
> end
> isitok=1;
> elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) &&
> (func(xnew(1,1),xnew(1,2))<values(sorted(1,1),:))
> disp('lets Contract with teta=beta...');
> xnew1=x(sorted(1,1),:)+(1+beta)*(xc-x(sorted(1,1),:))
> func(xnew1(1,1),xnew1(1,2))
> if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
> xnew=xnew1
> end isitok=1;
> end
> end
> x(sorted(1,1),:)=xnew;
> x
> disp('Finished================');
> x_archive=[x_archive;[x
> [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]]];
> end
>
> subplot(1,2,1)
> plot(x_archive(:,1),x_archive(:,2),'o-r')
> subplot(1,2,2)
> patch(x_archive(:,1),x_archive(:,2),x_archive(:,3),'r')
> view(3)
>
> values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
>
> sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
>
> disp(sprintf('Number of iteration= %10d\n',l))
> disp(sprintf('minimum of f(x) is %10.5f that occurs at
> x=[%10.5f,%10.5f]\n',min(values),x(find(values==min(values)),1),x(find(values==min(values)),2)))
>

Subject: optimization of fiber bragg grating using nelder-mead simplex

From: kuffal al

Date: 4 Apr, 2010 16:43:06

Message: 3 of 5

i actually dont really understand about the objective function bcoz there is a lots of equation. so which one should i choose or how should i create an objective function.


> You might want to use some built-in MATLAB functions.
> fminsearch performs the Nelder-Mead simplex algorithm.
> It accepts a single vector of inputs, say x = [x(1),x(2),x(3)], and
> works on your objective function that should take x and return f(x),
> where f(x) is the function to minimize. For your problem, take x =
> [L,A,dn], and write your objective function in terms of x.
>
> It seems to me that your code is doing some very inefficient
> root-finding as well. Take a look at the fzero function for a solid
> root-finding implementation.
>
> fminsearch reference:
> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fminsearch.html
>
> fminsearch example:
> http://www.mathworks.com/access/helpdesk/help/techdoc/math/bsgpq6p-10.html
>
> fzero reference:
> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fzero.html
>
> fzero example:
> http://www.mathworks.com/access/helpdesk/help/techdoc/math/bsgpq6q-44.html
>
> Alan Weiss
> MATLAB mathematical toolbox documentation
>
> kuffal al wrote:
> > can someone help me to optimize the fiber bragg grating(fbg) using
> > nelder mead simplex algoritm. i have both codes (both codes working
> > well) the fbg n nelder mead but i fail to integrate them. the parameter
> > that i want to optimize are the length (L), period (A) and the induced
> > index change (dn). so here are the codes
> >
> > (1)the fbg code
> > %==========================================================================
> > % matlab script for solving the coupled-mode equation using % the
> > transfer matrix method
> > %==========================================================================
> >
> > clear all
> > close all
> > clc
> >
> > %==========================================================================
> > % fibre simulation parameters
> > %==========================================================================
> >
> > L = 0.0100; % Grating Length
> >
> > dn = 1.0e-4; % Induced index change
> >
> > design_lambda = 1550e-9; % design wavelength
> >
> > %==========================================================================
> > % finding neff(effective reflective index) for the designed fibre
> > %==========================================================================
> >
> > format long
> > n=40;
> > a=0.5;
> > increment=2;
> > tolerance=1e-4;
> >
> > nco=1.447; % nco-core index
> > ncl=1.444; % ncl-cladding index
> > xstart=((ncl-0.0112)+nco)/2;
> > x=xstart;
> > dx=increment;
> >
> > %==========================================================================
> > % to find the root of neff(effective reflective index)
> > %==========================================================================
> >
> >
> > for m=1:n
> > neff= sqrt((ncl^2)+(a*x)*(nco^2-ncl^2));
> > while dx/x>tolerance
> > if neff~=sqrt((ncl^2)+(a*(x+dx))*(nco^2-ncl^2));
> > dx=dx/2;
> > else
> > x=x+dx;
> > end
> > end
> > route(m)=x;
> > dx=increment;
> > x=1.0001*x;
> > end
> > neff=route(40:end);
> >
> > V =(2*pi/design_lambda)*a*sqrt(nco^2-ncl^2); % normalize frequency
> >
> > A=design_lambda/(2*neff); % Grating period 5.36e-7;
> >
> > lambda=(1549.50e-9:0.005e-9:1550.50e-9); % Spectra range
> >
> >
> > %==========================================================================
> > % Coupled-mode Theory
> > %==========================================================================
> >
> > v=1; % visibility, v assumed to be 1
> >
> > k=(pi./lambda)*v*dn; % AC coupling coefficient
> >
> > % DC(period averaged) coupling coefficient
> >
> > % for stronger grating, sigma_dc assumed to be 0
> >
> > sigma_dc =0; %(2*pi./lambda).*dn;
> >
> > % general DC coupling coefficient
> >
> > sigma = 2*pi*neff.*(1./lambda - 1/(2*neff*A)) + sigma_dc;
> >
> >
> > %==========================================================================
> > % Power Reflection coefficients(R)
> > %==========================================================================
> >
> > K2 = k.*k;
> >
> > sigma2 = sigma.*sigma;
> >
> > difference = K2-sigma2;
> >
> > up = sinh(L*sqrt(difference)).^2;
> >
> > down = cosh(L*sqrt(difference)).^2 - sigma2./K2;
> >
> > R = up./down;
> >
> > Rt = transpose(R);
> >
> > %==========================================================================
> > % saving the points
> > %==========================================================================
> >
> > save 'fileToRead1.txt' Rt -ascii
> >
> > %==========================================================================
> > % plots for reflection spectra for Bragg reflector
> > %==========================================================================
> >
> > figure(1)
> >
> > plot(lambda,R,'b');
> >
> > axis( [1549.50e-9, 1550.50e-9, 0, 1]);
> >
> > title('Reflection Spectra of Bragg grating');
> >
> > xlabel('Wavelength(m)');
> >
> > ylabel('Reflectivity(p.u)');
> >
> > figure(2)
> > stairs(lambda,R,'r');figure(gcf)
> > axis( [1549.50e-9, 1550.50e-9, 0, 1]);
> >
> > --------------------------------------------------------------------------------------------------------------------
> >
> >
> > (2) the nelder-mead code
> >
> > %*******************************************
> > %* Nelder and Mead Technique *
> > %* Siamak Faridani *
> > %* University of Oklahoma *
> > %* 16 Oct 2005 *
> > %* www.PitchUp.com *
> > %*******************************************
> >
> >
> > clc;
> > clear all; close all;
> > N=2;
> > alpha=2;
> > beta=.25;
> > gamma=2.5;
> > Precision=.00001
> >
> > %Insert your starting point and the function here
> > x01=1.5505*1e-6;
> > x02=1.550e-6;
> > func=@Modified_fbg
> >
> >
> > %do not modify the following lines
> > %==============================
> >
> > x=[[x01,x02]];
> >
> >
> >
> > delta1=(((N+1)^.5+N-1)/(N*sqrt(2)))*alpha;
> >
> > delta2=(((N+1)^.5-1)/(N*sqrt(2)))*alpha
> > x=[x;[x(1,1)+delta1,x(1,2)+delta2];[x(1,1)+delta2,x(1,2)+delta1]];
> > x_archive=[x
> > [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]];
> > l=0
> > x
> > disp('start================');
> > Precision_inverse=Precision^-1;
> > while
> > ~(isequal(round(x(1,:)*Precision_inverse),round(x(2,:)*Precision_inverse),round(x(3,:)*Precision_inverse)))
> >
> > l=l+1;
> > values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
> >
> > sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
> >
> > xc=(1/N)*(x(sorted(1,2),:)+x(sorted(1,3),:));
> > disp('trying Normal Reflection');
> > xnew=x(sorted(1,1),:)+(1+alpha)*(xc-x(sorted(1,1),:))
> > isitok=0;
> > func(xnew(1,1),xnew(1,2))
> > while isitok==0
> > if ((func(xnew(1,1),xnew(1,2))>values(sorted(1,3),:)) &&
> > (func(xnew(1,1),xnew(1,2))<values(sorted(1,2),:)))
> > disp('f(l)<f(x_new)<f(g)')
> > isitok=1;
> > elseif (func(xnew(1,1),xnew(1,2))<values(sorted(1,3),:))
> > disp('lets Expand...');
> > xnew1=x(sorted(1,1),:)+(1+gamma)*(xc-x(sorted(1,1),:))
> > func(xnew1(1,1),xnew1(1,2))
> > if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
> > xnew=xnew1
> > end
> > isitok=1;
> > elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) &&
> > (func(xnew(1,1),xnew(1,2))>values(sorted(1,1),:))
> > disp('lets Contract with teta=-beta...');
> > xnew1=x(sorted(1,1),:)+(1-beta)*(xc-x(sorted(1,1),:))
> > func(xnew1(1,1),xnew1(1,2))
> > if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
> > xnew=xnew1
> > end
> > isitok=1;
> > elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) &&
> > (func(xnew(1,1),xnew(1,2))<values(sorted(1,1),:))
> > disp('lets Contract with teta=beta...');
> > xnew1=x(sorted(1,1),:)+(1+beta)*(xc-x(sorted(1,1),:))
> > func(xnew1(1,1),xnew1(1,2))
> > if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
> > xnew=xnew1
> > end isitok=1;
> > end
> > end
> > x(sorted(1,1),:)=xnew;
> > x
> > disp('Finished================');
> > x_archive=[x_archive;[x
> > [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]]];
> > end
> >
> > subplot(1,2,1)
> > plot(x_archive(:,1),x_archive(:,2),'o-r')
> > subplot(1,2,2)
> > patch(x_archive(:,1),x_archive(:,2),x_archive(:,3),'r')
> > view(3)
> >
> > values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
> >
> > sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
> >
> > disp(sprintf('Number of iteration= %10d\n',l))
> > disp(sprintf('minimum of f(x) is %10.5f that occurs at
> > x=[%10.5f,%10.5f]\n',min(values),x(find(values==min(values)),1),x(find(values==min(values)),2)))
> >

Subject: optimization of fiber bragg grating using nelder-mead simplex

From: James Allison

Date: 6 Apr, 2010 20:16:09

Message: 4 of 5

It seems that this might be formulated more naturally as a constrained
design optimization problem, i.e., you need to come up with both an
objective function AND constraint functions that fit the following form:

minimize f(x)
subject to g(x) <= 0
with respect to x

where x is the optimization vector (the independent variables you are
trying to adjust), and g may be a vector-valued function. What should be
the objective and what should be the constraints requires some
domain-specific knowledge about your application. What is most important
about the design? Are you trying to minimize or maximize something
(i.e., and objective)? Are there some limits that should not be exceeded
(i.e., constraints)?

fminsearch (Nelder-Mead) only applies to unconstrained problems, so it
won't work if you formulate your problem as a constrained optimization
problem. fmincon works well if the functions are smooth. patternsearch
would be the next thing to try if the functions are not smooth. The
genetic algorithm (ga) can solve very difficult optimization problems,
but requires a large number of function evaluations.

-James

kuffal al wrote:
> i actually dont really understand about the objective function bcoz
> there is a lots of equation. so which one should i choose or how should
> i create an objective function.
>
>
>> You might want to use some built-in MATLAB functions.
>> fminsearch performs the Nelder-Mead simplex algorithm.
>> It accepts a single vector of inputs, say x = [x(1),x(2),x(3)], and
>> works on your objective function that should take x and return f(x),
>> where f(x) is the function to minimize. For your problem, take x =
>> [L,A,dn], and write your objective function in terms of x.
>>
>> It seems to me that your code is doing some very inefficient
>> root-finding as well. Take a look at the fzero function for a solid
>> root-finding implementation.
>>
>> fminsearch reference:
>> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fminsearch.html
>>
>> fminsearch example:
>> http://www.mathworks.com/access/helpdesk/help/techdoc/math/bsgpq6p-10.html
>>
>>
>> fzero reference:
>> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fzero.html
>>
>> fzero example:
>> http://www.mathworks.com/access/helpdesk/help/techdoc/math/bsgpq6q-44.html
>>
>>
>> Alan Weiss
>> MATLAB mathematical toolbox documentation
>>
>> kuffal al wrote:
>> > can someone help me to optimize the fiber bragg grating(fbg) using >
>> nelder mead simplex algoritm. i have both codes (both codes working >
>> well) the fbg n nelder mead but i fail to integrate them. the
>> parameter > that i want to optimize are the length (L), period (A) and
>> the induced > index change (dn). so here are the codes
>> > > (1)the fbg code
>> >
>> %==========================================================================
>>
>> > % matlab script for solving the coupled-mode equation using % the >
>> transfer matrix method
>> >
>> %==========================================================================
>>
>> > > clear all
>> > close all
>> > clc
>> > >
>> %==========================================================================
>>
>> > % fibre simulation parameters
>> >
>> %==========================================================================
>>
>> > > L = 0.0100; % Grating Length
>> > > dn = 1.0e-4; % Induced index change
>> > > design_lambda = 1550e-9; % design wavelength
>> > >
>> %==========================================================================
>>
>> > % finding neff(effective reflective index) for the designed fibre
>> >
>> %==========================================================================
>>
>> > > format long
>> > n=40;
>> > a=0.5;
>> > increment=2;
>> > tolerance=1e-4;
>> > > nco=1.447; % nco-core index
>> > ncl=1.444; % ncl-cladding index
>> > xstart=((ncl-0.0112)+nco)/2;
>> > x=xstart;
>> > dx=increment;
>> > >
>> %==========================================================================
>>
>> > % to find the root of neff(effective reflective index)
>> >
>> %==========================================================================
>>
>> > > > for m=1:n
>> > neff= sqrt((ncl^2)+(a*x)*(nco^2-ncl^2));
>> > while dx/x>tolerance
>> > if neff~=sqrt((ncl^2)+(a*(x+dx))*(nco^2-ncl^2));
>> > dx=dx/2;
>> > else
>> > x=x+dx;
>> > end
>> > end
>> > route(m)=x;
>> > dx=increment;
>> > x=1.0001*x;
>> > end
>> > neff=route(40:end);
>> > > V =(2*pi/design_lambda)*a*sqrt(nco^2-ncl^2); % normalize frequency
>> > > A=design_lambda/(2*neff); % Grating period
>> 5.36e-7;
>> > > lambda=(1549.50e-9:0.005e-9:1550.50e-9); % Spectra range
>> > > >
>> %==========================================================================
>>
>> > % Coupled-mode Theory
>> >
>> %==========================================================================
>>
>> > > v=1; % visibility, v assumed to
>> be 1
>> > > k=(pi./lambda)*v*dn; % AC coupling coefficient
>> > > % DC(period averaged) coupling coefficient
>> > > % for stronger grating, sigma_dc assumed to be 0
>> > > sigma_dc =0; %(2*pi./lambda).*dn;
>> > > % general DC coupling coefficient
>> > > sigma = 2*pi*neff.*(1./lambda - 1/(2*neff*A)) + sigma_dc;
>> > > >
>> %==========================================================================
>>
>> > % Power Reflection coefficients(R)
>> >
>> %==========================================================================
>>
>> > > K2 = k.*k;
>> > > sigma2 = sigma.*sigma;
>> > > difference = K2-sigma2;
>> > > up = sinh(L*sqrt(difference)).^2;
>> > > down = cosh(L*sqrt(difference)).^2 - sigma2./K2;
>> > > R = up./down;
>> > > Rt = transpose(R);
>> > >
>> %==========================================================================
>>
>> > % saving the points
>> >
>> %==========================================================================
>>
>> > > save 'fileToRead1.txt' Rt -ascii
>> > >
>> %==========================================================================
>>
>> > % plots for reflection spectra for Bragg reflector
>> >
>> %==========================================================================
>>
>> > > figure(1)
>> > > plot(lambda,R,'b');
>> > > axis( [1549.50e-9, 1550.50e-9, 0, 1]);
>> > > title('Reflection Spectra of Bragg grating');
>> > > xlabel('Wavelength(m)');
>> > > ylabel('Reflectivity(p.u)');
>> > > figure(2)
>> > stairs(lambda,R,'r');figure(gcf)
>> > axis( [1549.50e-9, 1550.50e-9, 0, 1]);
>> > >
>> --------------------------------------------------------------------------------------------------------------------
>> > > > (2) the nelder-mead code
>> > > %*******************************************
>> > %* Nelder and Mead Technique *
>> > %* Siamak Faridani *
>> > %* University of Oklahoma *
>> > %* 16 Oct 2005 *
>> > %* www.PitchUp.com *
>> > %*******************************************
>> > > > clc;
>> > clear all; close all;
>> > N=2;
>> > alpha=2;
>> > beta=.25;
>> > gamma=2.5;
>> > Precision=.00001
>> > > %Insert your starting point and the function here
>> > x01=1.5505*1e-6;
>> > x02=1.550e-6;
>> > func=@Modified_fbg
>> > > > %do not modify the following lines
>> > %==============================
>> > > x=[[x01,x02]];
>> > > > > delta1=(((N+1)^.5+N-1)/(N*sqrt(2)))*alpha;
>> > > delta2=(((N+1)^.5-1)/(N*sqrt(2)))*alpha
>> > x=[x;[x(1,1)+delta1,x(1,2)+delta2];[x(1,1)+delta2,x(1,2)+delta1]];
>> > x_archive=[x >
>> [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]];
>> > l=0
>> > x
>> > disp('start================');
>> > Precision_inverse=Precision^-1;
>> > while >
>> ~(isequal(round(x(1,:)*Precision_inverse),round(x(2,:)*Precision_inverse),round(x(3,:)*Precision_inverse)))
>> > > l=l+1;
>> >
>> values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
>> > >
>> sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
>> > > xc=(1/N)*(x(sorted(1,2),:)+x(sorted(1,3),:));
>> > disp('trying Normal Reflection');
>> > xnew=x(sorted(1,1),:)+(1+alpha)*(xc-x(sorted(1,1),:))
>> > isitok=0;
>> > func(xnew(1,1),xnew(1,2))
>> > while isitok==0
>> > if ((func(xnew(1,1),xnew(1,2))>values(sorted(1,3),:)) && >
>> (func(xnew(1,1),xnew(1,2))<values(sorted(1,2),:)))
>> > disp('f(l)<f(x_new)<f(g)')
>> > isitok=1;
>> > elseif (func(xnew(1,1),xnew(1,2))<values(sorted(1,3),:))
>> > disp('lets Expand...');
>> > xnew1=x(sorted(1,1),:)+(1+gamma)*(xc-x(sorted(1,1),:))
>> > func(xnew1(1,1),xnew1(1,2))
>> > if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
>> > xnew=xnew1
>> > end
>> > isitok=1;
>> > elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) && >
>> (func(xnew(1,1),xnew(1,2))>values(sorted(1,1),:))
>> > disp('lets Contract with teta=-beta...');
>> > xnew1=x(sorted(1,1),:)+(1-beta)*(xc-x(sorted(1,1),:))
>> > func(xnew1(1,1),xnew1(1,2))
>> > if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
>> > xnew=xnew1
>> > end
>> > isitok=1;
>> > elseif (func(xnew(1,1),xnew(1,2))>values(sorted(1,2),:)) && >
>> (func(xnew(1,1),xnew(1,2))<values(sorted(1,1),:))
>> > disp('lets Contract with teta=beta...');
>> > xnew1=x(sorted(1,1),:)+(1+beta)*(xc-x(sorted(1,1),:))
>> > func(xnew1(1,1),xnew1(1,2))
>> > if (func(xnew1(1,1),xnew1(1,2))<func(xnew(1,1),xnew(1,2)))
>> > xnew=xnew1
>> > end isitok=1;
>> > end
>> > end
>> > x(sorted(1,1),:)=xnew;
>> > x
>> > disp('Finished================');
>> > x_archive=[x_archive;[x >
>> [func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))]]];
>> > end
>> > > subplot(1,2,1)
>> > plot(x_archive(:,1),x_archive(:,2),'o-r')
>> > subplot(1,2,2)
>> > patch(x_archive(:,1),x_archive(:,2),x_archive(:,3),'r')
>> > view(3)
>> > > values=[func(x(1,1),x(1,2));func(x(2,1),x(2,2));func(x(3,1),x(3,2))];
>> > >
>> sorted=[find(values==max(values)),6-find(values==max(values))-find(values==min(values)),find(values==min(values))];
>> > > disp(sprintf('Number of iteration= %10d\n',l))
>> > disp(sprintf('minimum of f(x) is %10.5f that occurs at >
>> x=[%10.5f,%10.5f]\n',min(values),x(find(values==min(values)),1),x(find(values==min(values)),2)))
>> >

Subject: optimization of fiber bragg grating using nelder-mead simplex algoritm

From: shammosa@gmail.com

Date: 5 Apr, 2013 22:22:08

Message: 5 of 5



calculation first...please respond

please help, i hope you can send me the reference of effective idex equatio that you had used in the calculation first...please respond

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us