// Slower C++ solution:
#include <bits/stdc++.h>

#define N 100005

#define ll long long

#define For(i,x,y) for(int i=(x);i<=(y);++i)

#define Rof(i,x,y) for(int i=(x);i>=(y);--i)

#define p_b push_back

#define pii pair<int,int>

#define mp make_pair

#define fr first

#define sd second

#define mod 1000000007

using namespace std;



ll _2[N],_2_[N],fac[N],invf[N];

int nxt[N],id[N],len[N],deg[N],a[N],tot[N],cir[N],d[N],pos;

bool vis[N],inq[N];

vector<pii > v[N];



void dfs(int x,int dep){

	if(inq[x] || cir[x]){

		if(!cir[x]){

			len[x]=dep-d[x],pos=x,cir[x]=x,id[x]=1;

			if(d[x]-1) v[x].push_back(mp(1,d[x]-1));

		} else if(dep-1) v[cir[x]].push_back(mp(id[x],dep-1));

		return; 

	}

	d[x]=dep,inq[x]=1,vis[x]=1;

	dfs(nxt[x],dep+1),inq[x]=0;

	if(pos && dep>d[pos]) cir[x]=pos,id[x]=id[nxt[x]]+1;

}

int ___(int x,int y){ if(x>y) return 2; if(x<y) return 0;return 1; }

ll qpow(ll x,int y){

	ll tmp=1;

	for(;y;y>>=1,(x*=x)%=mod) if(y&1) (tmp*=x)%=mod;

	return tmp;

}

int main(){

	int n;ll ans=1,ans2=1;

	scanf("%d",&n);

	_2[0]=_2_[0]=1;ll fuck=qpow(2ll,mod-2);

	fac[0]=fac[1]=invf[0]=invf[1]=1,_2[1]=2,_2_[1]=fuck;

	For(i,1,n) scanf("%d",&a[i]),nxt[i]=a[i],deg[nxt[i]]++;

	For(i,2,n){

		_2[i]=(_2[i-1]+_2[i-1])%mod;

		_2_[i]=_2_[i-1]*fuck%mod;

		fac[i]=1ll*fac[i-1]*i%mod;

		invf[i]=1ll*(mod-mod/i)*invf[mod%i]%mod; 

	}

	For(i,1,n) (invf[i]*=invf[i-1])%=mod;

	For(i,1,n)

		if(!deg[i]) pos=0,dfs(i,1);

	For(i,1,n) if(!vis[i]) pos=0,dfs(i,1);

	For(i,1,n){

		if(deg[i]>2) return puts("0"),0;

		else if(!cir[i] && deg[i]>=2) return puts("0"),0;

	}

	For(i,1,n) sort(v[i].begin(),v[i].end());

	For(i,1,n){

		if(!v[i].size()){ if(len[i]) tot[len[i]]++;continue; }

		if(v[i].size()==1){

			if(len[i]<v[i][0].sd) return puts("0"),0;

			else if(len[i]>v[i][0].sd)(ans*=2ll)%=mod; 

		} else{

			For(j,0,(int)v[i].size()-2){

				ll p=1ll*___(v[i][j+1].fr-v[i][j].fr,v[i][j].sd);

				if(!p) return puts("0"),0;

				else (ans*=p)%=mod;

			}

			int p=v[i].size()-1;

			if(len[i]-v[i][p].fr+v[i][0].fr<v[i][p].sd) return puts("0"),0;

			else if(len[i]-v[i][p].fr+v[i][0].fr>v[i][p].sd) (ans*=2ll)%=mod;

		}

	}

	For(i,1,n){

		if(!tot[i]) continue;

		ll qwq=0;

		For(j,0,tot[i]/2){

			int x=j*2;

			ll tmp=qpow(1ll*i,j)*fac[tot[i]]%mod*invf[tot[i]-x]%mod*_2_[j]%mod*invf[j]%mod;

			(tmp*=_2[((i%2) && i!=1)*(tot[i]-x)])%=mod;

			(qwq+=tmp)%=mod;

		}

		(ans2*=qwq)%=mod;

	}

	printf("%lld\n",ans*ans2%mod);

}

// Faster C++ solution for the same problem (measured speedup: 3.586961835820978x):
#include <bits/stdc++.h>

using namespace std;



typedef long long ll;

const int maxn = 100005, mod = 1e9 + 7;



inline int gi()

{

	char c = getchar();

	while (c < '0' || c > '9') c = getchar();

	int sum = 0;

	while ('0' <= c && c <= '9') sum = sum * 10 + c - 48, c = getchar();

	return sum;

}



inline int add(int x, int y) {return x + y >= mod ? x + y - mod : x + y;}



int n, ans, a[maxn], cir[maxn], flen[maxn], sum[maxn], f[maxn], deg[maxn], vis[maxn];



void solvecir(int x)

{

	int len = 0, s = 0, lst = 0, slen = 0;

	while (cir[x]) {

		++len; cir[x] = 0;

		if (flen[x]) {

			if (!s) s = lst = len, slen = flen[x];

			else {

				ans = (ll)ans * ((flen[x] < len - lst) + (flen[x] <= len - lst)) % mod;

				lst = len;

			}

		}

		x = a[x];

	}

	if (!s) ++sum[len];

	else ans = (ll)ans * ((slen < len - lst + s) + (slen <= len - lst + s)) % mod;

}



void solve()

{

	for (int i = 1; i <= n; ++i) {

		if (deg[i]) continue;

		int x = i, len = 0;

		while (!cir[x]) x = a[x], ++len;

		flen[x] = len;

	}

	ans = 1;

	for (int i = 1; i <= n; ++i) if (cir[i]) solvecir(i);

	for (int i = 1; i <= n; ++i) {

		if (!sum[i]) continue;

		f[0] = 1;

		for (int j = 1; j <= sum[i]; ++j) {

			if (i > 1 && (i & 1)) f[j] = add(f[j - 1], f[j - 1]);

			else f[j] = f[j - 1];

			if (j > 1) f[j] = (f[j] + (ll)f[j - 2] * (j - 1) % mod * i) % mod;

		}

		ans = (ll)ans * f[sum[i]] % mod;

	}

}



int main()

{

	n = gi();

	for (int i = 1; i <= n; ++i) a[i] = gi(), ++deg[a[i]];

	for (int i = 1; i <= n; ++i) {

		if (vis[i]) continue;

		int x = i; while (!vis[x]) vis[x] = i, x = a[x];

		if (vis[x] != i) continue;

		while (!cir[x]) cir[x] = 1, x = a[x];

	}



	for (int i = 1; i <= n; ++i)

		if ((cir[i] && deg[i] > 2) || (!cir[i] && deg[i] > 1)) return puts("0"), 0;



	solve();

	printf("%d\n", ans);



	return 0;

}