Matlab数值解下面的方程,具体见图：如有建设性思路、想法，或具体的程序一定加分我这样这样画图的：fun_fan.m文件中：function [ fun ] = fun_fan( d,N_eff )lambda_0=1.55e-6;k_0=2.*pi./lambda_0;epsilon_m=-127.8-3.

Matlab数值解下面的方程,具体见图：如有建设性思路、想法，或具体的程序一定加分我这样这样画图的：fun_fan.m文件中：function [ fun ] = fun_fan( d,N_eff )lambda_0=1.55e-6;k_0=2.*pi./lambda_0;epsilon_m=-127.8-3.
Matlab数值解下面的方程,具体见图：

如有建设性思路、想法，或具体的程序一定加分

我这样这样画图的：

fun_fan.m文件中：

function [ fun ] = fun_fan( d,N_eff )

lambda_0=1.55e-6;

k_0=2.*pi./lambda_0;

epsilon_m=-127.8-3.21i;

epsilon_d=1;

alpha_xd=(N_eff.*k_0.^2-epsilon_d.*k_0.^2).^0.5;

alpha_xm=(N_eff.*k_0.^2-epsilon_m.*k_0.^2).^0.5;

fun=epsilon_m.*alpha_xd.*tanh(1/2.*alpha_xd.*d)+epsilon_d.*alpha_xm;

end

test.m文件中：

clear all;

clc;

ezplot(@(d,N_eff)fun_fan(d,N_eff),[0,400e-9,1,2]);%限定范围

%ezplot(@(d,N_eff)fun_fan(d,N_eff));%不限范围

title('N_{eff} vs. d');

这样不能得到Neff实部与d的关系图，Neff虚部与d的关系图，请大牛指导，万分感谢。。。

Matlab数值解下面的方程,具体见图：如有建设性思路、想法，或具体的程序一定加分我这样这样画图的：fun_fan.m文件中：function [ fun ] = fun_fan( d,N_eff )lambda_0=1.55e-6;k_0=2.*pi./lambda_0;epsilon_m=-127.8-3.
我觉得首先要定义neff的实部和虚部
function [ fun ] = fun_fan( d,N_real,N_image)
N_eff=N_real+i*N_image;
lambda_0=1.55e-6;
k_0=2.*pi./lambda_0;
epsilon_m=-127.8-3.21i;
epsilon_d=1;
alpha_xd=(N_eff.*k_0.^2-epsilon_d.*k_0.^2).^0.5;
alpha_xm=(N_eff.*k_0.^2-epsilon_m.*k_0.^2).^0.5;
fun=epsilon_m.*alpha_xd.*tanh(1/2.*alpha_xd.*d)+epsilon_d.*alpha_xm;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all;
clc;
ezplot(@(d,N_real)fun_fan(d,N_real,0),[0,400e-9,1,2]);%限定范围
figure
ezplot(@(d,N_image)fun_fan(d,1,N_image))
%ezplot(@(d,N_eff)fun_fan(d,N_eff));%不限范围
title('N_{eff} vs.d');
貌似对于给定的d,n_real和n_image没有唯一解,所以可能在画neff_real to d的图时,要更改iamge的值.

难啊

我编程解了一下那个方程组，并画出了neff的实部和d的关系图，发现d在0~400e-9之间变化时，neff的值基本不变，因为你的d的变化范围实在太窄了，在这么窄的范围内解基本稳定！下面是我变得程序：
huatu.m文件中：
x0=[5,5,5];
X=ones(9,3);
options=optimset('MaxIter',5e6,'TolX',1e-9);

全部展开

我编程解了一下那个方程组，并画出了neff的实部和d的关系图，发现d在0~400e-9之间变化时，neff的值基本不变，因为你的d的变化范围实在太窄了，在这么窄的范围内解基本稳定！下面是我变得程序：
huatu.m文件中：
x0=[5,5,5];
X=ones(9,3);
options=optimset('MaxIter',5e6,'TolX',1e-9);
for i=1:9
d=50e-9*(i-1);
X(i,:)=fsolve(@(x)myfun1(d,x),x0,options);
end
neff=X(:,1)
plot(0:50e-9:400e-9,real(neff))
%%%%%%下面是myfun1.m函数文件
function F=myfun1(d,x)
neff=x(1);
alphaxd=x(2);
alphaxm=x(3);
lambda0=1.55e-6;
k0=2*pi/lambda0;
epsilonm=-127.8-3.21i;
epsilond=1;
F=[alphaxd^2-neff*k0^2+epsilond*k0^2;alphaxm^2-neff*k0^2+epsilonm*k0^2;
epsilonm*alphaxd*tanh(0.5*alphaxd*d)+epsilond*alphaxm];

收起